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Abstract: The dynamic renormalization group (RG) is used to study the large-distance and long-time limits of viscous 
and resistive incompressible magnetohydrodynamics subject to random forces and currents. The scale-dependent viscosity and 
magnetic resistivity are derived and used for carrying out RG-improved perturbation theory. This is applied to derive both 
the asymptotic scaling and the overall proportionality coefficients for both the velocity and magnetic field correlation functions 
as well as the kinetic and magnetic energy density spectral functions. The Kolmogorov, Iroshnikov-Kraichnan, as well as 
other energy spectra, formally can be obtained by suitable choice of injected noise, although the method limits the validity of 
these energy spectra only to the asymptotic regime. Injection of a random magnetic helicity is considered, its RG-improved 
spectral density derived, and its contribution to the velocity and magnetic field correlation functions determined. The RG 
scaling solutions are used to determine information at asymptotic scales about energy and helicity cascade directions and 
mixing between magnetic and kinetic energy. Some of the results found here also are shown to be valid for the Navier-Stokes 
hydrodynamic equation. 

PACS: 95.30.Qd, 05.10.Cc, 47.27.Gs, 98.62.En 

I. INTRODUCTION 

The nonlinear equations of magnetohydrodynamics (MHD) admit turbulent solutions that can cover a wide range 
of spatial scales. One expects to see the turbulent regime of MHD as the most probable state of a magnetofluid 
present in physical systems ranging from the solar corona (Van Ballegooijen 1985; 1986) and solar wind (Burlaga 
1991) to galactic magnetic fields (DeYoung 1980; 1992). MHD turbulence in the early universe may have even left a 
detectable present-day signature at the cosmological scale (Brandenburg et al. 1996). 

To study questions such as the asymptotic scaling behavior in magnetic field correlations, the presence or absence 
of kinetic and magnetic energy cascades, the helicity spectrum and its corresponding cascades, as well as the effect 
of random perturbations, one needs to consider the evolution of the magnetofluid over a wide range of spatial and 
temporal scales. In general, turbulence is an essential feature in coming to terms with these problems, whose rigorous 
treatment calls for an analysis of fully nonlinear MHD. 

For the regime that interests us, the large-distance and long time behavior of the fluid velocity and magnetic fields, 
the equations of MHD lend themselves to an analytical treatment by applying the methods provided by the dynamic 
renormalization group (RG) . The aim of the RG is to describe quantitatively how the dynamics of a system evolves 
as one changes the length and time scales of the phenomena under study. The RG is an excellent tool for exploring 
the physics of complicated physical systems that are characterized as having many interacting degrees of freedom 
and are subject to fluctuations covering many different length scales. RG methods were first used in the study of 
hydrodynamic (non-magnetic) turbulence by Forster, Nelson and Stephen (1977). They studied the large-distance 
and long-time properties of velocity correlations generated by the Navier Stokes equations for a randomly stirred 
incompressible fluid with a variety of Gaussian noise spectra. More recently these methods were used by Yakhot 
and Orszag (1986) and Dannevik, Yakhot and Orszag (1987) in an attempt to predict several parameters in fully 
developed hydrodynamic turbulence. A detailed critique of their use of the RG can be found in the paper by Eyink 

(1994) . Some recent reviews of the RG in hydrodynamic turbulence may be found in Frisch (1995) and McComb 

(1995) . 

There have been rather fewer applications of the RG to magnetized fluids, probably due to the complexity of the 
problem and the lack of experimental data for turbulent magnetic fields in laboratory plasmas. Indeed, as indicated 
above, MHD turbulence is much more likely to be of interest at astrophysical and cosmological scales. The previous 
works, that apply RG to MHD turbulence for addressing various specific problems, have been carried out by Fournier, 
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Sulem and Pouquet (1982), Longcope and Sudan (1991) and more recently by Camargo and Tasso (1992). A recent 
survey of MHD turbulence can be found in the monograph by Biskamp (1993). 

In this paper, RG will be applied to examine the asymptotic properties of turbulent MHD with a stochastic force 
term. Our analysis treats the full MHD equations, and thus differs from Longcope and Sudan (1991), who treated 
a type of "reduced MHD". Also, we treat all nonlinear terms as given in the MHD equations on an equal footing, 
in contrast to Fournier et. al. (1982), who weighted the inertial nonlinearity and Lorentz force differently. In this 
respect, the starting equations for us are the same as those of Camargo and Tasso (1992), although our treatment 
differs significantly from their work. In regards to the calculation technique, we will cast the stochastic MHD equations 
into a functional path integral formalism, in contrast to the direct iterative solution of the equations of motion carried 
out in (Camargo and Tasso 1992). The functional path integral formalism permits use of useful field theory methods, 
which are much more compact and efficient than direct iterative expansions. Although in this paper we use this 
formalism only to implement our one-loop RG analysis, the formalism can be applied much more generally to study 
stochastic MHD. With respect to results, we differ from all the previous works (Fournier et al. 1982; Longcope and 
Sudan 1991; Camargo and Tasso 1992). We will examine cascade directions and energy mixing, present solutions of 
the correlation functions in the scaling regime and implement a RG improved perturbation theory. All of these results 
are being explored and presented here for the first time. Furthermore, we also will implement a RG analysis when 
magnetic helicity is present, and to our knowledge no such purely analytic treatment of this problem has been carried 
out before in the literature. 

The paper is organized as follows. In Section || we introduce the stochastic equations of incompressible MHD. These 
equations then are written in terms of Elsasser variables, which make manifest the symmetry between the velocity 
and magnetic field for incompressible MHD (Elsasser 1950). This symmetry is useful in checking the correctness of 
the renormalization group equations that will be derived and for understanding the nature of the correlation functions 
and energy spectra for the velocity and magnetic fields. In the stochastic equations, the stochastic random forces 
and currents are intended to model both random initial conditions and subsequent disturbances and fluctuations. The 
noise therefore serves a two-fold purpose and it summarizes at a coarse-grained level our ignorance about the details 
of the numerous small-scale phenomena. For example, as applied to the problem of galactic or even cosmic fields, 
the noise would account for the many diverse particle physics mechanisms (e.g., early phase transitions, vacuum 
bubble collisions, fluctuating Higgs gradients, superconducting cosmic strings, etc.) that have been proposed for the 
generation of seed magnetic fields (for a review see eg. Enqvist (1998)). In this particular context then, the noise 
simulates in part the effects of the particle-physics degrees of freedom at the micro-scale. Although we will not 
examine models to motivate the noise forces, we emphasize that the noise need not be a stirring force (in the strict 
sense of turbulence modelling) but also can account, in part, for real physical fluctuations and perturbations intrinsic 
to the system. 

In Section [II we turn to the renormalization group analysis of incompressible MHD. First the physical motivations 
are discussed for implementing the RG techniques to MHD. Then Subsec. [II A examines the naive scaling properties 
of stochastic MHD under independent space and time rescalings. Here, two fundamental scaling exponents are 
intro duced, that encode all the information pertaining to the asymptotic behavior of the solutions of MHD. Subsec. 
[II B describes the two step RG procedure of coarse graining and rescaling. It then explains the relation between scale 
invariant solutions of the MHD equations and fixed points of the RG analysis. In particular, at a fixed point the MHD 
equations become scale invariant, and a unique solution exists for the pair of scaling exponents, which are calculable 
by the RG procedure. To fix ideas, a simple example of the RG procedure is given in the context of linearized MHD. 
Of course this exampl e doe s not account for the nonlinear interactions, so crucial for establishing and maintaining 
turbulence. In Subsec [II C the RG analysis is done for the full MHD equations, in particular including the nonlinear 
terms. As alluded to above, most of the applications of RG to hydrodynamics and MHD require one to iterate the 
equations of motion to a particular order in the nonlinear coupling. This can result in rather cumbersome expressions 
after a few iterations. Here, we prefer to exploit the streamlined functional integral and allied diagrammatic methods 
which allow for the rapid identification of the one-particle irreducible (1PI) diagrams whose renormalization then yields 
directly the independent RG equations for the parameters appearing in the MHD equations. Since these techniques 
lie somewhat out of the mainstream of the paper and can be applied more generally to the problem, we relegate 
them to a series of technical Appendices, where the complete renorma lizati on of the response function, noise spectral 
function and interaction vertex is carried out in full detail. In Subsec III C mainly the differential RG flow equations 
are presented, and the calculations leading to these equations are outlined with references to the appropriate places 
in the Appendices for the various details. 



In Sec. IV, the RG flow equations are analyzed with a pair of dimensionless couplings. Then in this section, the 



fixed points are calculated, their stability properties are elucidated and the associated critical exponents are obtained 
for each fixed point. All these results are given for any spatial dimension d and for a general class of stochastic noise 
forces, characterized by a single spectral force exponent y, that ranges from short (y < 0) to long (y > 0) range. 
In Section ^ we derive and solve differential equations for the scale dependent viscosity and resistivity, derive the 
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scaling forms obeyed by the velocity and magnetic field correlation functions, and derive both the scaling form and 
RG improved form of the kinetic and magnetic energy spectra. Here, an expression is said to be renormalization group 
improved if the bare parameters in that expression are replaced by their scale-dependent forms, typically calculated 
to some order in perturbation theory. Thus, an improved quantity, whether it be a correlation function, an energy 
spec trum , etc., is one which combines perturbation theory with the renormalization group (Hochberg et al. 1999a). In 



Sec. VI, we turn to consider the injection of random magnetic helicity, which is modeled with a power law spectrum, 
and compute the RG-improved helicity spectrum. In this Section the correlation functions are recomputed, now taking 
into account the helicity contribution. Then once again, these correlation functions are computed in the RG improved 



form. An alternative way to understand scaling and the approach to scaling is given in Section VII, which ma kes u se 
of the Callan-Symanzik equation whose solution yields the improved correlation function automatically. Sec. VIII is 
concerned with energy and helicity cascades in hydrodynamics and MHD. A methodology for extracting information 
about cascades from our RG calculations is developed and then implemented. The predictions from our RG analysis 



are given for cascades and cascade directions at the asymptotically largest length scales. Finally, in Se ction IX we 
draw conclusions and summarize our work. Here some discussion is provided, extending from Subsect VC, on the 
mixing of magnetic and kinetic energy at asymptotic scales. 

As mentioned above, the Appendices contain most of the technical calculations of our RG analysis as well as a 
derivation of the functional path integral formalism that is generally applicable to stochastic MHD. In particular, 
the Appendices are organized as follows. In Appendix A we set up the dynamic functional for incompressible non- 
relativistic stochastic MHD and use this to derive the basic elements needed for a Feynman diagram representation 
of the solutions of the MHD equations and define the response function, noise spectral function and vertex function. 
In Appendix B we carry out the momentum-shell integration needed to renormalize the response function. From this 
we obtain the one-loop corrections to the kinematic and magnetic viscosities. The one-loop renormalization of the 
noise spectral function is carried out in Appendix C. For y > — 2, which corresponds to noise spectra from long range 
down to reasonably short range, we find that there is in fact no one-loop correction to the noise amplitude in the 
hydrodynamic limit. The one-loop renormalization of the vertex function vanishes identically for all gaussian forcing 
functions that are temporally white, but for arbitrary spatial correlations, as demonstrated in Appendix D. Some 
useful identities needed for averaging products of unit vectors over the unit (i-sphere are collected in Appendix E. 

II. STOCHASTIC MHD EQUATIONS 

The equations of motion for viscous and resistive, incompressible, non-relativistic, magnetohydrodynamics are given 
by (Biskamp 1993) 

r)v 1 1 

— + (vV)v = — V(p + — ) + — (B-V)B + ^V 2 v + ff v , (1) 
at p v 87T ' inp 

and 

^(B'V)v-(vV)B + £v 2 B + I / B , (2) 



subject to the constraints 



and 



V-v = 0, (3) 



V • B = 0, (4) 

where v denotes the fluid velocity, B the magnetic field, p the fluid pressure and p the density. The speed of light is 
c and a is the electrical conductivity of the magnetofluid. We allow for the presence of random forces and random 
currents, denoted by ff v (x,t) and and f]B(x,t), respectively. We assume these forces are Gaussian with zero mean, 
(ffv) = {^b) = 0. The random vector forces and currents are specified by their two-point correlation functions, which 
for convenience we specify directly in Fourier space and in terms of their individual vector components as 

((&)„(£, w)(^) m (fcV)> = 2D v g(k)(2TT) d+1 V nm {k)5{Lu + uj')8 d {k + k'), (5) 
((f)B)n(k,u)(ffB) m (k',u')) = 2D B g(k)(27T) d+1 P nm (k)S(Lu + Lu')S d (k + k'), (6) 
((T&)„(£,w)(fo) m (£V)>=(). ( 7 ) 
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Here, d is the number of space dimensions, the indices are n,m — 1, 2, • • • , d, and the angular brackets denote the 
average taken over the fluctuations. The presence of both kinetic and magnetic helicities can be taken into account 
and modeled by means of suitably defined noise terms. We defer the treatment of random helicity to a later section. 
The transverse projection operator P m „(fc) is needed to ensure the noise is solenoidal, i.e., k ■ ff v = k ■ ffs = 0, and 
hence compatable with incompressibility (V • v = 0) and the absence of magnetic monopoles (V • B = 0). The Fourier 
transform of the random forces is given by 

ffv,B(x,t) = J ^ J ^f) ViB (k,u) exp(ik-x-iut), (8) 

and this is the convention used throughout this paper. The magnitudes of the random forces and currents are given 
by D v > and Db > 0, respectively, and so provide a measure of the size of the fluctuations. The power law function 
g(k) — k~ y with exponent y is intended to describe the noise spectrum in the inertial range. Eq. (Q) states that 
the random currents and forces are uncorrelated. The noise exponent y will be treated as a model parameter. Note 
depending on the sign of y, the noise fluctuations are either more correlated at short distances (y < 0) or at large 
distances (y > 0). 

We now cast the above coupled vector stochastic partial differential equations into a manifestly symmetric form 
first by defining B — * H/y/Anp, p — > p/p and = (p + \B 2 ). Then, following Elsasser (1950), we add and subtract 
the two equations of motion Eqs. ([[]) and @j to obtain 

dP - - - 

— + (Q ■ V)P = -Vp* + 7+ V 2 P + 7 -V 2 Q + ff P , (9) 



at 

^ + (P ■ V)Q = -Vp, + 7 +V 2 Q + 7 -V 2 P + f) Q , (10) 



where 



P = v + B, Q = v-B, (11) 

1+ = \{v + vb), 7- = ^(v-Vb), (12) 

rjp = rfv +Vb, ff Q =ff v - tJb- (13) 

In terms of these Elsasser variables, the noise statistics is completely determined from Eqs. (||) - (Q) together with 
the definitions Eq. dl3|) as 



{{rf P ) n (k,u){vp) m (k' =2[D v + D B }g(k)(2TT) d+1 P nm (k)S(u; + cu / )d d (k + k') 7 

= ((rf Q ) n (k,u J )(rf Q ) m (k',u J ')) 1 (14) 
((Vp)n(k, io)(rf Q ) m (k', a/)) = 2[D V - D B ] g(k)(2ir) d+1 P nm (k)6(co + u')8 d (k + k 1 ). (15) 

We henceforth write vb = ■r £ —, to indicate the magnetic resistivity. This pair of dynamical equations Eqs. (||) and 
( |l0| ) is symmetric under the interchange P <-> Q and ffp <-> rjQ- Furthermore, in the limit of vanishing magnetic field, 
P = Q and vanishing magnetic noise, ffp — Vq, the pair reduces to two identical copies of the Navier Stokes equation 
with a random noise source. These symmetries and limits are useful in checking results derived from these equations 
and can be used as a diagnostic in comparing our calculations with other RG results obtained from MHD and from 
incompressible hydrodynamics. 

For incompressible MHD, we can take the divergence of both equations Eqs. ([|) and (|l0|) to eliminate the pressure 
term Vp* since V ■ v = V ■ B = 0. This yields 

p* = -^d n d k {QkPn), (16) 

which can be regarded as the equation of state for this system. We note that in solving for , there are no boundary 
terms that contribute to Eq. ([l6|). As in Eqs. (0) and (|J), the stochastic random forces and currents can be assumed 
to be solenoidal, since any longitudinal components can be absorbed into the pressure . This simplifies the equations 
even further. To do this, we introduce the transverse projection operator, which is the same one appearing above in 
the noise correlation functions 



n = (Sjn - d j^2 d n), (17) 
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for j,n = 1, 2, • • • , d. Of course, this projector is a nonlocal operator in coordinate space (as it involves the Green 
function V -2 ), but its Fourier transform is local and easier to handle in actual computations. Transforming Eq. (]lT 
using Eq. M), gives 



Pjn(fc) — ( Sj n 



kj k n 
k 2 

Now can be eliminated from Eqs. (^|) and (|lo| ) to give, in coordinate space, 

BP 



BQj 
Bt 



X P jn di(PiQn) - 7+V 2 Q, - 7-V 2 P, = (vq), 



(18) 

(19) 
(20) 



This is the final form of the stochastic MHD equations that serves as the starting point for our renormalization group 
analysis. The presence of the projection operator makes it preferable to express these equations in component form, 
as we have done here. We introduce a parameter Aq = 1 in front of the nonlinear terms, which is useful for organizing 
perturbative expansions of the solutions of these equations. This is merely a bookkeeping device, since our true (and 
dimensionless) perturbation expansion parameter will be defined and calculated below. Finally, note that correlations 
(and other functions) in the fluid velocity and magnetic fields can be computed in terms of these Elsasser fields, since 
the two sets of variables are linearly related as 



v(x,t) = ^[P(x,t) + Q(x,t) 
B(x,t) = ~ [P(x,t) - Q(x,t) 



(21) 
(22) 



III. RENORMALIZATION GROUP (RG) ANALYSIS OF MHD 



The renormalization group (RG) is a systematic method for studying a system at different length scales. One of the 
key motivations for implementing this method, and the main application in this paper, is to verify a scaling hypothesis. 
In vague terms a scaling hypothesis is related to the statement that in the asymptotic regime, the correlations in 
a system all can be expressed in terms of a single suitable length scale. For the renormalization group method, 
vindication of a scaling hypothesis is related to showing that the underlying dynamical equations are self-similar at 
different scales, up to an overall rescaling of the space-time coordinates and the dynamical variables of the system. 
Sometimes self-similar solutions are referred to as fixed points, since the parameters of the system do not change 
under the RG transformation. For hydrodynamics and magnetohydrodynamics, the scaling hypothesis that has been 
considered (Forster et al. 1977; Yakhot and Orszag 1986; Danncvik et al. 1987; Fournier et al. 1982; Longcope and 
Sudan 1991; Camargo and Tasso 1992) and will be considered here, is associated with the asymptotic regime of large 
length and time separations in the correlation functions. 

In this Section, the above statements will be explained in detail and then applied. Subsection [II A first will 



explain the scale transformations of the MHD equations that are important for the renormalization gro up pr ocedure. 
Subsection III B then explains the renormalization group in the context of MHD. Finally Subsection [II C presents 
the results of applying the renormalization group to MHD. This last subsection requires the technically most involved 
aspect of the calculation in this paper. For ease of presentation, all the technical steps have been given in the 
Appendices and Subsec. Ill C simply quotes the results and gives some guidance to the Appendices. 



A. Scale Transformation 



In preparation for a dynamical renormalization group analysis, and as an important and necessary preliminary step, 
the stochastic MHD equations Eqs. (|l^) and ( pp| ) are submitted to independent global rescalings of space and time 
as 

x = sx' 1 => V=-V', 

s 

- m = '-'m> (23) 
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and the (Elsasser) fields as 



P{x, t) = 8*P{?, t'), Q{x, t) = s x Q'(f , t'), (24) 

where s > 1. To fix ideas, one can think of this rescaling as examining the original unprimed system in a new primed 
coordinate system, with the two related by the above rescalings. 

The above transformation involves two exponents, the dynamic exponent z and the "roughness" exponent x- 
The choice of the above scale transformation is made because in the renormalization group analysis, the number 
of independent space and time rescalings one must consider depends on the degree of symmetry and anisotropy of 
the system under study. As we are considering homogeneous and isotropic magnetohydrodynamic turbulence, this 
scaling introduces two a-priori independent exponents. Alternatively, for example in reduced (2+l)-dimensional MHD 
where one has a dominant background magnetic field, say along the z-direction, one scales the (x, y)-plane and the 
perpendicular z directions independently (Longcope and Sudan 1991) and this introduces three independent scaling 
exponents: z, Xi £• 

An important step in examining self-similar behavior in the renormalization group approach is to examine the 
original unprimed equations in the primed coordinate system and compare the form of the equations. The precise 
relevance of this step will be explained in the next subsection. However, at a procedural level, in performing this step, 
one term in both set of equations always can be made similar by dividing through by an appropriate factor. For us, 
we will choose this term to be the time derivative terms dP/dt 7 dQ/dt. Thus Eqs. ( |l9| ) and (20) when expressed in 
the primed coordinates are 

d tl P\(x' ,t') + s z+ x- 1 \ P' jn d' l (Q' l Pn) = s z - 2 7 +V' 2 Pj(f ,*') - s z - 2 j-V 2 Q' j (x' ,0 

+ a *-x+&v- d -') T iP i (5f,1?), (25) 



and similarly for the Q-equation. There are three points about Eq. ( |25| ) to note. First, in arriving at Eq. ( p5|) a factor 
of s~ z+x was divided through the equation so that the time derivative terms in this and Eqs. (|19| ) and ( p0| ) are the 
same. Second, the transformation behavior of the noise functions is fully specified by the space and time rescalings 
in Eq. ( p3| ) and the form of the noise correlations in Eqs. (||) - ((?]). Third, the MHD equations generically have an 
implicit short distance cut-off Ax > Axo or in momentum space |k| < A ~ 1/Axo. For distances less than this cut-off 
the hydrodynamic approximation is invalid in properly representing the physics. The continuum description of the 
fluid breaks down. The scale transformations above also act on this short distance cut-off so that in the primed frame 
this cut-off becomes Ax' — sAxo or, equivalently A' = A/s. Since s > 1, we see that the cut-off (or lattice spacing) 
has been scaled up (in real space) towards the infrared. 

We compare this scaled equation of motion Eq. ( p5| ) with its unsealed version Eq. (|l9|). If the scaled equation 
Eq. (|2^) is to describe the same dynamics, albeit at the larger scale (s > 1), then this implies that the coefficients 
appearing in it must change accordingly. By "same dynamics" we mean of course that the mathematical form of 
the dynamical equation does not change as one changes the scale of the problem. This implies the following scaling 
behavior of the parameters appearing in the MHD equations: 

A — S Ao, 

7^=* 2 - 2 7-, 
A' = s z ~ 2x ~ d+y A 

B 1 = s z - 2x - d+y B. (26) 

We have defined the noise-amplitude combinations A = D v + Db and B = D v — Db (Thus, the zero magnetic noise 
limit Db — is had by taking A = B). 



It is important to appreciate that Eq. ( |25| ) with the definitions Eq. (|26|) and Eqs. (19) and ( |2C| ) are simply 
two expressions of the same equation. In particular no requirement of scale invariance is made. We simply have 
required that any changes they suffer due to this scaling transformation be absorbable into redefinitions of the (bare) 
parameters appearing in the equations. 

There is one final notational comment in regards to the scale transformation. Due to the central role that this 
transformation plays in the renormalization group analysis, it is standard to assign all quantities in the dynamical 
equation a scaling dimension. Typically the momentum coordinates, which scale as s _1 , is referred to as having 
scaling dimension 1, so that length has scaling dimension — 1. Similarly time, the fields, and all the parameters in Eq. 
( pri| ) have scaling dimension given by minus the exponent to which s is raised. 
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B. Renormalization Group Procedure 



The RG transformation is now rather standard fare and is discussed in a number of excellent texts and monographs 
(Ma 1976; Amit 1978; Zinn- Justin 1996; Cardy 1996). Here we briefly review the RG procedure. The RG is a 
symmetry transformation that carries a system from one length scale to another through a two step procedure, coarse 
graining and rescaling. As with any symmetry transformation, such as rotation, translation etc., it is not necessary 
that the system under examination be invariant to this transformation. However, should this be the case, then is as 
also the case with any symmetry transformation, important information about the system can be deduced without 
detailed knowledge of its entire dynamics. In particular, invariance under the RG transformation will imply scale 
invariance in the system. The existence of such an invariance implies the behavior of correlation functions at different 
scales can be understood relative to their behavior at one given scale. 

Let us review the two step RG procedure. Consider the exact MHD Eqs. (Q) and ( |To| ) or equivalently Eqs. ( |Tgj ) 
and ( po| ) . As stated earlier, there is an implicit assumption that these equations are only valid down to some distance 
scale x > Axo or in momentum space k < A ~ 1/Axq. This is the smallest length scale at which the MHD dynamics, 
Eqs. (^) and ([lO]), are meant to be valid. Given this starting point, the RG is implemented in the following two steps: 

Coarse-graining: The first step is to solve the dynamics at the very shortest length scales. In Fourier space we 
integrate the solutions of the stochastic equations of motion over a thin momentum shell A/s < |fc| < A where 
s = 1 + S and < 6 << 1. Physically, this step serves to integrate out the fast or short wavelength modes and 
thins out the degrees of freedom. This elimination process leads to "running" in the parameters (which depend on 
,s = e l ) appearing in the MHD equations. The short distance cut-off or lattice spacing gives rise to a corresponding 
momentum cut-off A/s. The resulting MHD equations after this step will be given by 

d t p< + \<(e)p rn d l (Q<p<) - i${£)v 2 p< 

- -/< {t)V 2 Qf - (r}^)j + less relevant = 0, (27) 

and likewise for Qj. Clearly in performing this step additional terms will appear in the dynamic equations. These 
terms have not been explicitly written above, but have been referred to suggestively as "less relevant", since we will 
see later that for the leading asymptotic regime they will not be important. Nevertheless, at this stage if we retain 
all the terms in the coarse grained equations Eq. (^7|), including the "less relevant" terms, the dynamics from these 
equations is exactly the same as from the original MHD Eqs. (||) and (|Io|), for computing correlations for length scales 
larger than Ax > 2irs/A. In particular, correlations computed from Eqs. @ and ( fl0| ) or Eqs. ( p7j ) will be the same, 
i.e. 

(P i (x 1 ,t 1 )P j (x2,t2)) = (P^{x 1 ,t 1 )P 3 <(x2 > t 2 )), etc. (28) 

Rescaling: The second step is to represent the above coarse grained Eqs. (^7|) in the primed coordinates based on 
the rescalings in Eqs. ( p3f ) and ( p^ ) with s = e l . In the primed coordinates, the MHD equations become 

d t ,p< + x'MPjndKQW) - y + (£)V 2 p; 

- 7^(^)V' 2 Qj- - s {y - z - d),2 {rip) J {x' , t') + less relevant rescaled = 0, (29) 
and likewise for Q. In these equations 





= ^ +x " 1 A <, 






7-W 


z-2 < 

= s 7_ , 


A'(£) 


— gZ-2-X-d+y 


B'(£) 


— s z-2x-d+y q< 



(30) 

Note, if we solve Eqs. ( |29| ) and examine correlations at a distance d — (x' x — x' 2 ), this is equivalent to solving for the 
same correlation in the original Eqs. (||) and (JTo|) , at distance e e d = (x\ — x 2 )- 

At this point we can consider Eqs. (^9|) in place of the original MHD Eqs. (||) and (0), and repeat these two 
steps iteratively. A very interesting thing can happen upon sufficient iteration, which is a key practical importance of 
the RG approach. The set of parameters in iterative versions of Eqs. ( p9| ) could be the same after a certain number 
of successive RG transformations. These points in the parameter space are called fixed points. At the fixed point, 
correlations computed from either of two successive equations clearly are identically the same, since the successive 
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equations are exactly the same. This fact has nontrivial consequences, since due to the scaling relations between the 
two sets of equations, it also means we can learn about the behavior of correlations at different scales of the original 
MHD equations. In particular consider the field P( x, t), P'(x',t') (likewise Q) from two successive iterations of Eqs. 
(^9|). If these two successive versions of Eqs. (|29| ) are exactly the same then it is an identity that (dropping the 
indices) 

(P(iMi)P(*2,ta)) = (P(^iA)P(^t' 2 )), (31) 

and likewise for Q. However, based on the rescaling Eqs. ( p3| ) and (|24|) we also have 

(P(x 1 ,t 1 )P(x 2 ,t 2 )) = e 2 ^{P{e- l x u e- zt t l )P{e- i x 2 , e^ta)), (32) 

which therefore gives a relation between correlations at two different length/time scales. 

To help focus these ideas, let us consider a simple example of a fixed point solution in RG for the case of MHD 
without the nonlinear terms. In this case, the MHD dynamics from Eqs. (^]) and (|lfj| ) are 

dtPj - 7+V 2 P, - 7-V 2 Q, - = 0, (33) 



d t Qj ~ 7+V 2 Q, - 7 -V 2 P, - (vq)j = 0. (34) 

Implementing the RG procedure, the first step is coarse graining. Since this is a free (and hence, linear) theory, 
this step is trivial. We simply remove the highest modes between A/s < |fc| < A from the P and Q fields. The second 
step then is rescaling which based on Eqs. ( p3| ) leads to 

d v p;{x', o - a *- 2 7+ v ,2 p; - s z - 2 7 _v ,2 q;. - s *-*+^- d -*V. ( x ',t') = o, (35) 

and similarily for Q. Observe if 

z = 2 

y-d + z 

X = 2 ' ( 36 > 

then both Eqs. ( |33] ) and ([55]) are the same, and so this is a fixed point of the RG. 

Using this fact to compute correlations at two different scales in the momentum space correlation functions yields 

(P(fci s wi)P(A& s W2)) = e 2(x+d+2) VV^e z V)PV^,e^ 2 )) 

= e^ +d +^(P(e%, e*V)P(e%, e^ 2 )>, (37) 

and similarily for Q. Since this is a free theory, we can check this result of the RG procedure from the exact correlation 
functions. For this note that the above fixed point Eq. (|3^) is the same for both linearized MHD Eqs. ( |33| ) and ([34] ) 
or linearized hydrodynamics, which arises from these equations by identifying P = Q. Also, the exact solutions of 
both linearized MHD and linearized hydrodynamics have the same basic form, but the latter are much less elaborate. 
Thus, it is sufficient for our demonstrative purpose here to quote the simpler case of linearized hydrodynamics, for 
which the exact two-point correlation function is 



, ,->,, , k 1 y S d (k 1 +k 2 )d(u 1 +uj 2 ) 
(P(fci,wi)P(fc 2 ,w 2 )) ~ const x -J 2 , / U 2\2 ' ( 38 ) 



where v = j + +7- from Eq. (12). The scaling transformations Eqs. (23) and (pjj)) applied to this imply 

^(^(Pfe^^e'VlPfe^.e^)) ~ e 2{ x+ d +z -i- y -d- z ) mnst x 

k^ v 5 d (k 1 +k 2 )5(ui +lu 2 ) 



2\2 



(39) 



Comparing the scaling behavior from the RG analysis Eq. fl37|) and the exact solution, Eq (|39|), we see they agree. 

In summary, returning to the opening statements of this section, finding fixed points in the RG procedure corre- 
sponds to identifying scale invariance in the MHD equations with respect to the scale transformations Eqs. ( p3|) and 
(p4|). As one practical matter, since the iterative procedure in Eqs. (^) and ( |29[ ) is focussed on the behavior of the 
parameters, it is standard and convenient to express the RG procedure through a set of differential equations for the 
evolution of the parameters, eg. dX(£)/d£, dj+(£) / dl, etc. The fixed points then correspond to the solutions where all 
parameters have zero evolution, that is, d\(£)/d£ = dj + (£)/d£ = ■ ■ ■ = 0. 
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C. Differential Renormalization Group Equations for MHD 



This subsection presents the RG evolution equations for MHD. There are several steps to this calculation, all of 
which have been carried out in full detail in the Appendices. Here these steps only will be outlined with reference to 
the appropriate Appendices for the complete details, and then the results will be given. For solving the equations of 
stochastic MHD, in Appendix A they have been Fourier transformed and cast into a dynamic functional formalism. 
This formalism allows a systematic loop-expansion, which in other words is a perturbation expansion in the amplitude 
of the random fluctuations. The basic building blocks for constructing the perturbative expansion are expressed 
diagrammatically through the Feynman rules in Fig. 1. Based on this approach, the information contained in the 
explicit solutions of P, Q is equivalently contained, but much more compactly, in the response function, noise spectrum 
and vertex function, whose definitions and diagrammatic representations (expressed to one-loop order) are given in 
Figs. (2a)-(2c). Thus, we work directly with these objects as they will yield, to a given order in perturbation theory, 
the corrections to the MHD parameters Ao,7+,7-,^4 and B. As stated in the last subsection, these corrections 
result from moving to ever larger length scales while taking into account the fluctuations from the smaller scales. 
We have calculated the one loop corrections to the response function, the noise spectral function and the nonlinear 
vertex function in the Appendices B, C and D respectively, and from them, have obtained the corresponding one- loop 
corrections to the parameters Ao, 7+, 7-, A and B in the hydrodynamic limit (i.e., the limit of large distance and long 
time). We point out that in arriving at the RG flow equations below, we have taken the limit lo — ► at the outset 
of the calculations as described in the Appendices. This is equivalent to taking the limit t — > 00, and so serves to 
"project" out the asymptotic time limit. Finally the second step of the RG procedure is to rescale the coarse grained 
quantities based on Eqs. (p3h and (p~~ 



Applying these steps to the one- loop expressions Eqs. (D8), (C9), (CK), ( B17 ), and (B18), computed in the 
Appendices, we obtain the following set of one-loop renormalization group equations for nonrelativistic incompressible 
MHD: 

d\ 
dl 
dA 

~dJ 
dB_ 

H 



(*+ 


x- 


1)A , 


(*- 




-d + y)A, 


(*- 


2x- 


-d + y)B, 


(7+ 


+ 7- 


){(z-2) 


(7+ 


- 7- 


){(*-2) 



dl w+ '-' \ y ' ' 2 u a (7++7-) 3 

- h + - {(. - 2) + l^-'^L} , (40, 



where 



_ (d*-y-4)S d 

Ad ~ d(d + 2)(2n)d ' (41) 

and Sd — Y(d/2) • Here, S — 4 is the numerical symmetry factor associated to the diagram representing the one- loop 
correction to the response function (see, Fig. (2a)), which has been calculated in Appendix B. The short-distance 
or ultraviolet cut-off (or, lattice spacing) is A, below which the hydrodynamic or continuum description of the fluid 
breaks down. We c an s et A = 1 without loss of generality. In fact, none of the RG results depend explicitly on this 



cut-off. In Section VII we provide an alternative way to understand scaling which makes implicit use of the short- 
distance divergences that result from formally taking the limit A — > 00 in the one-loop expressions. It is clear that Eq. 
( fi"o| ) describe the scale dependence of the parameters of MHD. Their solution therefore yields the effective parameters 
Ap (i) , A (£) , B (I) , 7 + (I) and J-(i)- Substitution of these scale-dependent parameters into the MHD equations Eqs. 
( |l9| ) and (|2(j) yields the effective dynamical equations corresponding to this scale as parametrized by I > 0. 
In terms of the dimensionless couplings defined by 



9± - |^A^- 4 ( {A ^ B \ 3 , (42) 



the RG flow can be summarized by means of the two equations 







-j^=g + ([A-d + y]-3Sg + ), (43) 

^=g-([i-d + y]-3Sg.), (44) 

which shows that flow takes place in a two-dimensional coupling or parameter space with axes (g_|_,<?_). Here, 5 = 4 
is the same numerical symmetry factor as mentioned above. 

By shutting off the one- loop corrections (by putting Ao = 0) in these equations Eq. ((40|), one immediately recovers 
the naive scaling in Eq. (p6|). This demonstrates the importance of the nonlinear interaction for scaling behavior. We 
see that the scale dependence of the viscosity and resistivity is modified from its naive behavior. However, the noise 
amplitudes and the formal expansion parameter do not get modified at one-loop; see Appendices C and D for further 
details. 



IV. FIXED POINT SOLUTIONS AND EXPONENTS 

The solutions of the differential RG equations Eq. (f4(i|) (or equivalently Eqs. ( f43| ) and (|44|)) tell us how the 
parameters appearing in the MHD equations Eqs. ( [l9| ) and ( p0[ ) change with length and time scale. At a fixed point, 
as the name implies, the parameters no longer change under a further RG transformation, and one is instead in the 
scaling regime, where the dynamics is scale invariant. To find the renormalization group fixed points, we set the 
differential RG equations to zero, which yields algebraic equations. 

The RG equations Eqs. p3] ) and (^) have the two fixed point solutions: 

(t) (g* + ,gl) = (0,0) and, (45) 

(5 +' 5 - )= ( i2 ' 12 J" (46) 

We wish to emphasize the fact that g± provides the bona- fide expansion parameter: it is dimensionless and is less 
than one whenever (4 — d + y) < 12. This latter condition is satisfied for all the noise spectra and spatial dimensions 
considered in this paper. The stability properties of each fixed point is determined by performing a linear analysis 
about each point. That is, we expand the couplings g± = g± + Sg± where Sg± is a small fluctuation, and then 
substitute this expansion into Eqs. (^3|) and (^4|) to obtain the RG equations for the linearized perturbations about 
the fixed points Eqs. @ and @. This yields 

(<) = (4 - d + y)Sg + + 0{5g\), (47) 

(ii)^j± = -{±-d + y)8g + + 0{8gl). (48) 

Identical equations hold for Sg^ by simply replacing 8g + — > Sg- in these equations. These results tell us that when 
(4 — d + y) > 0, the trivial fixed point Eq. ( [45| ) is infrared unstable and repulsive (since the perturbation grows for 
£ — ► oo) while the non-trivial fixed point Eq. (^6|) is infrared stable and attractive. What this means physically is 
the following. We imagine probing the system at a given initial length and time scale. There is a corresponding 
set of parameters Ao, A, B, 7 + , 7^ associated to this scale. Geometrically, this can be represented by a point in (a 
five-dimensional) parameter space. Then, tracing the effective dynamics of the system at ever larger length and time 
scales generates a directed path in this space which starts off from this initial point. This flow therefore is always 
driven away from the repulsive fixed point and towards the attractive fixed point, which therefore acts as a "sink". 
Therefore, the asymptotic dynamics is uniquely determined by the attractive fixed point. If (4 — d + y) < 0, then the 
stability properties of these two fixed points are interchanged. In this case, all RG flows, irrespective of where they 
start off in parameter space, inevitably end up in the neighborhood of the trivial fixed point Ao = 0. This implies 
the asymptotic dynamics is linear (and trivial). We see the stability properties of fixed points depend on the spatial 
dimension d and on the noise spectrum through the noise exponent y. For the Kolmogorov spectrum y — d and for the 
Iroshnikov-Kraichnan spectrum y = d — 1/4, the trivial fixed point Eq. ( |45| ) is always repulsive while the non-trivial 
fixed point Eq. ( f46|) is always attractive in the infrared. Similar conclusions hold for a white noise spectrum, y = 0, 
provided the dimension d < 4. Finally, the fixed point critical exponents are obtained by substituting the values of 
the fixed points into the original set of RG equations Eq. ( |4C| ) , which yields the solutions 

(z,x)= (2,~(2 + 2/-d)) , for fixed point (i), (49) 

(z, X )= (~(2 + d- y),~(l+y -dj) , for fixed point (ii). (50) 
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The exponents for the trivial fixed point are simply the generalization of those of the Edwards- Wilkinson model of 
surface growth subject to spatially correlated noise (Barabasi and Stanley 1995). Note furthermore that the exponents 
satisfy \ + z — 1 a t the nontrivial fixed point (i.e., when Ao ^ 0). This in fact had better be the case as dictated by 
the renormalization group equation for Ao in Eq. (|40| ) which admits a nontrivial fixed point solution iff z + x — 1 = 0. 
So, this serves as an important consistency check on our calculation of the exponents. Finally, note that for either 
fixed point, the dynamic exponent is non-zero, z =/= 0. 



V. SCALING AND THE APPROACH TO ASYMPTOTIA 



With the fixed point solutions obtained in the previous sections we now turn to the asymptotic scaling behavior 
(in the infrared) of the velocity and magnetic field correlation functions and energy density spectra. 



A. Renormalized scale dependent viscosity and resistivity 

The approach to the (nontrivial) fixed point can be quantified and is very useful for carrying out RG- improvement. 
The effective viscosity v(k) and magnetic resistivity ^s(fc) at the scale k result from all turbulent motions with wave 
numbers above k (i.e., at shorter length scales). To derive the renormalized viscosity and resistivity, we integrate Eqs. 
( pl7| ) and p!8| ) over a finite band of momenta from the cut-off A to A/s, but without rescaling. 



From taking first the sum and then the difference of these expressions, we derive the following infinitesimal recursion 
relations 

dv = 2A d (A + B) f A ^ (<*-»-*) . 



dt v 2 \e e J 

dv B _ 2A d (A - B) / A\ (d-v-i) 



(52) 



where we have set Ao = 1 and have used s — e . We assume that the noise amplitudes, A, B have achieved their 
asymptotic (i.e., near the fixed point) values. 



Integrating Eq. (51) yields the relation 



- 3 W = - 3 (°) + ^y-Q {(--^) d - y - 4 A^- 4 ) , (53) 

where ^(0) is the initial value of the viscosity. When a large range of wavenumbers is eliminated this takes a scaling 
form independent of initial properties and the ultraviolet cut-off A, 



( fc ) = ( 6 4+ A tf )V (W)/3 , (54) 



where k = e e A and A + B = D v > is the magnitude of the random force fluctuations. 

By a similar argument, we find that the k — ► scaling form of the magnetic resistivity is given by 

where A — B = Db > is the magnitude of the random current fluctuations. 

Thus we see that the viscosity and magnetic resistivity tend to increase as k — > (that is, for large length scales) 
provided (4 + y — d) > 0, which is the condition that the RG trajectories flow towards the nontrivial fixed point. 
This trend of increasing effective viscosity has also been reported as a prediction of the RG in pure Navier-Stokes 
(Dannevik et al. 1987) and in incompressible MHD (Longcope and Sudan 1991; Camargo and Tasso 1992), where in 
addition there is an increasing effective magnetic resistivity. These scale dependent quantities will be used repeatedly 
later on when we carry out the renormalization group improvement of the correlation functions and spectral densities. 
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B. Scaling form of the correlation functions 



The scaling properties of the P- field correlations, the Q-field correlations as well as the cross correlations provide 
direct information on the large-distance and long-time behavior of the fluid velocity and magnetic field via the simple 
linear rel ation s writt en in Eqs. ( pl| ) and (S3) . These correlations form part of a larger block- matrix structure as given 
in Eqs. (A26), (A29) and (A27), but we can easily work out the scaling properties of the individual elements of this 



block array and further note from direct inspection of those expressions that (PjP n ) = (QjQ n ) and (PjQ n ) = (QjPn)- 
So, there arc really only two independent types of field correlation functions to consider, whether these are expressed 
in terms of the Elsasser variables or in terms the velocity and magnetic fields. In this Section, we derive the exact 
asymptotic scaling form of the correlation functions up to (as of yet) an unknown scaling function. This illustrates the 
power and limitations of pure scaling arguments. After the inclusion of magnetic hclicity, we return to the correlation 
functions and will determine their scaling functions explicitly by means of RG-improved perturbation theory. 

It will prove extremely useful to have the scaling form of these correlations available in both the space and time 
domain as well as in the Fourier (momentum and frequency) domain. In the scaling regime, the P-field scales as 



Pj(x, t) = s x P'^x", t') = s x P^s- 1 x, s- z t) = s x Pjis^x, s- z t), 



(56) 



and similarly for Q. The middle expressions are direct consequences of the scale transformations Eqs. ([23|) and ( |24| ) 
whereas in the scaling regime, the MHD equations are self-similar and then the rightmost equality also holds. Thus 
in the scaling regime the P — P correlation function scales as 



(Pj(£,t)P n (0,0)) = s 2x (P J -(s- 1 f ! s- z t)P„ (0,0)), 
C PjPn (x,t) = s*xCp jPn {8- x 2,8-*t). 



(57) 



Here by C PjPn , we mean the correlation function between the jth and the nth components of the P-field. Note 
that because of space and time translational invariance, the correlation for the fields depend only on the spatial 
separation x and temporal separation t. The above equality relates the correlation function at two different scales, 
whose separation is measured by s > 1. We emphasize that Eq. ([57]) is exact, when the system is in the scaling 
regime, that is to say, in the vicinity of a fixed point. 

We now will see what further information can be extracted from the expression Eq. (|57|). There is an implicit 
tensor structure in these functions, due to the vector indices, which we make explicit below in the Fourier domain. 
Consider spatial correlations at equal times (i.e., t = 0) then from Eq. (p7|) we have that 



C P]Pn (x,0) = \x\ 2 *C PjPn {l,Q), 
whereas temporal correlations at the same spatial point (i.e., x = 0) scale as 

Cp jPn (0,t)=t 2 ?C PjPn (0,1), 



(58) 



(59) 



where we have chosen s = \x\ in the former and s — t 1 ^ in the latter. Here, C P Pn (1, 0), Cp p„ (0, 1) are just constant 
tensor coefficients. The general scaling form of the correlation function valid for arbitrary spatial and temporal 
arguments can be written as 



C P]Pn (x,t) = \x\ 2x C PjPn 



t 



(60) 



where the scaling function C obeys the following limits: 



u -> 
u — > oo 



C PjPn {u) ^ C PjPn (1,0), 
Cp jP Ju)^C PPn (0,l)u- 



(61) 
(62) 



Note that Eqs. ( p8|) and ( |59| ) are just special cases of Eq. ( |60[ ). 

Identical results hold for the Q — Q correlation function by simply replacing the index labels P — > Q from Eq. (|57 
to Eq. (|62[). For the cross-correlations P — Q, just replace one of the P symbols in Eq. (|57| ) to Eq. ( |6^ ) by a Q. 
This is as far as one can go employing scaling arguments. However, this is enough for determining how the power-law 
correlations depend on the exponents z and which are obtained from solving for the RG fixed points. 

Fourier transforming Eq. ( |5^ ) yields equivalent information, but expressed instead in the momentum and frequency 
domain 
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CPjPnik^) = k 

with of course, the scaling functions Cpp = Cqq. 
independent scaling function, and thus we have 



-ax-* p. 



(63) 



The scaling of the cross-correlation function involves an a-priori 



C 



PjQ 



Xk,oj) = k- d -^- z v ]n {k)c PQ [- 



(64) 



Note that Eqs. (|63j) and (|64J) show that the tensor structure of the correlation functions is carried by the projection 
operator, as was to be expected. 

It should be clear from these considerations that the asymptotic behavior of the correlation functions and the 
quantities derived from them, arc all completely determined in terms of the two exponents z and x calculated above. 
Thus, using Eqs. ( pl| ) and (22) to transform back to the physical velocity and magnetic fields, we immediately find 
in the space and time domain 



t 



C 



P 3 Q n 



t 



(Vj(x,t)v n (0,0)) = |x| 2 * |C Pj p,. v ( 
{ Vj (x,t)B n (0,0)) = 0. 



(65) 

(66) 
(67) 



Note there is no net averaged cross-helicity. This is consistent with the fact that no helicity has been injected into 
the system. 

From Eqs. ( p35"| ) - (|67]), the spatial dependence of both fluid velocity and magnetic field correlations in the vicinity 
of the non-trivial IR stable fixed point scale as 



(v 3 (x,0)v n (0,0)) ~ (B 3 (x,0)B n (0,0)) ~ ri^"^, 



(68) 



where r = \x\ and the roughness exponent x = ^(1 + J/ — d) at the infrared stable fixed point Eq. ([46|). Note that for 
the Kolmogorov spectrum, y = d, and we recover the two-thirds law (r%) for the velocity correlation function in fully 
developed hydrodynamic turbulence (Frisch 1995) (and we predict that the spatial dependence of the magnetic field 
correlation function is a power-law obeying a two-thirds law for this noise spectrum). For the Iroshnikov-Kraichnan 
spectrum, the spatial correlation scales instead as r 3 . For the temporal correlations, from Eqs. ( |65|) - (|67j ) they scale 
as 



(« i (o,i)«„(o,o)) ~ (^-(o^^^oD-^^iw 



(69) 



where the dynamic exponent z = g(2 + d — y). For the Kolmogorov spectrum, the temporal correlations scale as ~ t, 
whereas for the I-K spectrum, they scale as ~ i 6 / 7 . By marked contrast, a white noise spectrum y = yields spatial 
correlations that scale as r~z and temporal correlations that scale as t~ "S in d — 3 dimensions. 

To help understand the above results about the correlation functions in the scaling regime, the case of d = 3, 
y = — 2 hydrodynamics will be computed explicitly. The above scaling relations of the correlation functions all apply 
to hydrodynamics by setting the magnetic field to zero, B = 0. For this case 4 — d + y < 0, thus the trivial fixed point 
is stable and so it is simple to obtain the exact correlation functions by inverse Fourier transforming Eq. (pq) to yield 



Cpp ^ (x, t) — const, x 



CXp( -4^ } 



/{vt) 



3/2 



const, x 



(M!)3/a exp(-J^) 
y vt ' K Avt 1 



(70) 



Let us verify this solution is consistent with the general scaling forms given above. Since at this trivial fixed point 
for d — 3, y = — 2 from Eq. (|49| ) z = 2, x = — 3/2, it follows that the above solution is consistent with the general 
scaling forms Eqs. ( p7| ) and (|60[). Furthermore, from the general scaling form Eq. (59) one expects Cpp{Q, t) ~ i~ 3 / 2 
which is consistent with the above solution. For t = 0, the inverse Fourier transform is ill-defined, so comparison is 
not relevant. This is not inconsistent with the general scaling form Eq. d58| ) since these general considerations do not 
guarantee that the function Cpp(l, 0), for example in Eq. (pq), necessarily is well defined. 
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C. Energy spectrum 



A quantity of fundamental importance in turbulence research is the energy spectrum which contains, in part, 
information regarding the flow of energy across different length scales. The energy spectrum also is key in defining 
and calculating the spectral energy density, which is the energy density per mode or wavenumber. The transfer of 
energy from smaller to larger (or larger to smaller) length scales as a system evolves in time is known as an inverse 
(direct) energy cascade. We can investigate whether and under what conditions, randomly forced stochastic MHD 
in the absence of magnetic helicity, subject to the noise spectrum of Eqs. (||) and (||), will exhibit energy cascades, 
and more generally, what is the asymptotic behavior of the energy spectrum. To do so, we will begin by focusing 
attention on the quadratic energy densities for the Elsasser fields and make direct use of the asymptotic fixed point 
scaling properties for the correlation functions deduced in the previous subsection. After this is done, we then turn to 
a complete calculation of the spectral energy density functions for both the physic al flu id velocity and magnetic field 



by means of improved perturbation theory using the results obtained in Subsection V A for the renormalized viscosity 
and magnetic resistivity. This will lead to identical scaling relations as obtained from the "pure" or naive scaling 
arguments, but with the added advantage of yielding the explicit form of the individual scaling functions, valid to 
the one-loop order in perturbation theory at which we are working. Moreover, the example will serve to demonstrate 
how RG-improvement works in practice. 

Consider therefore the total energy Epp in the Elsasser P field averaged over the random fluctuations (this will 
serve as a useful guide when we come to treat all the fields together below) 



E PP = \j d d x(P(x,tyP(x,t)), 
= itr j d d x{P J (x,t)P k (x,t)), 



= ~tr| d d x j J ^ d eMi{k + k')-x){P 3 {k,t)P k (k\t)), 

-^J^- d (P^t)P k (-k,t)), 

1 f d d k f duj f duj' , , .. . , , , -> ... 
= 2 tT J 2^ J ^ eM-K" + "')t)(Pi(k,u)P k (-k,u% 

= Ztr [** [ d ^Cp Pk (k^), (71) 

2 J {2-n) d J 2tt 3 kK ' h V ' 

where we have used (Pj(k,ui)P k (k' ,u>')) = (2ir) d+1 5 d (k + k')6(uj+uj')Cp j p k (k,uj), and the d-dimensional spatial volume 
factor is V = (2ir) d S d (6). From the identity in the last line of Eq. @, we can immediately identify the energy 
density (Epp/V), which is given by the indicated double integral of the correlation function. The spectral energy 
density, i.e., the energy density per unit wavenumber, is therefore given by 

1 f°° 
Epp(k) - [-^/(27r) d+1 ] k^ 1 / dto tr C P .p k (k,cj), (72) 

* J— oo 

since clearly 

/>oo 

Epp = V dkEpp(k). (73) 
J o 

By making use of Eq. (|63|) which holds in the scaling regime, we can discover how the spectral energy density scales 



as a function of the injected random noise spectrum. Inserting the expression Eq. (63) into Eq. (f72j) yields 



(d-l) S d 
2(2tt)° 



3C 



E PP( k ) = k ° I duC P p{u). (74) 



This follows after defining the exponent 



a = -l-2 X =-~ + ~(d-2/), (75) 
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and making use of the exponent identity z + \ = 1 , which holds at the non-trivial infrared stable fixed point Eq. (|4^ 
The steps leading to the derivation of the spectral energy density function EQQ(k) for the Q-field are identical to those 
outlined above, and in fact, we have Epp(k) — EgQ(k), as a consequence of the equality between the scaling functions 
Cpp = Cqq. We remark at this point that the choice of noise exponent y — d formally leads to the Kolmogorov 
energy density spectrum, since this yields a = — |, and the spectral energy density then scales as 

Ep P {k)^k'i. (76) 

By way of contrast, at the unstable and trivial fixed point, the spectral density would scale instead as 

E PP {k) ~ k~ 3 (77) 

for y — d. The difference in the scaling is due to the presence Eq. ((76|) or absence Eq. (^) of the nonlinear terms in 
the equations of motion which act to mix distinct modes. 

We note of course that the overall magnitude of the spectral density Eq. ( f74|) is determined in part by the scaling 
function Cpp. It has been common practice in many, if not all, RG analyses of turbulence to deduce only the 
power law scaling exponent of the physical quantities of interest while leaving the overall proportionality factors 
undetermined. This could lead to the impression that the RG is incapable of determining these overall factors. In 
fact, the renormalization group not only yields a determination of the scaling exponents but also provides quantitative 
information regarding the overall constant factors. This is achieved to a given order in perturbation theory by means 
of RG-improved perturbation theory (Hochberg et al. 1999a). We will make use of this technique to complete our 
discussion of the spectral energy density. We now turn to this task, which will be carried out in two steps. First we 
calculate the energy corresponding to the free- field limit, where the correlation functions are known explicitly, and 
then we use RG to improve the resulting energy by use of our perturbation theory. This results in the corrected form 
of the spectral energy, to the same order of perturbation theory at which we are working. We begin again at the 
level of the Elsasser fields and then at the end of the calculation make use of the linear transformation back to the 
physical fluid velocity and magnetic fields. The net result will be a one-loop expression for both the fluid velocity and 
magnetic field spectral energy densities. The calculation of Epp(k) carried out in detail above will serve as a useful 
template for this purpose. 

We define an energy function involving both the Elsasser fields P and Q at once (which is actually a 2 x 2 matrix, 
as we indicate with the dyad) 

- V f d d k f°° duj .» , r ,. 
E = TT tr / / 7T [C (k,L>)] mn , 



2 J (2ir) d J_ 00 2n 

poo 

V dkE(k), (78) 



where the energy spectral density function (also a 2 x 2 matrix) is given by 

E (k) = [-S d /(27r) d+1 ]k d ~ 1 / dw tr [C (k, u)] mn . (79) 

This expression is exact. We can calculate it to one-loop order by means of RG-improved perturbation theory. This 
proceeds by the two above mentioned steps. That is to say, first we evaluate this spectral density using the free- field 
correlation function (denoted with the zero subscript) 

[Co (k {k)go{k,u)D {k)g (-k,-u), (80) 

and then improve the resulting expression by replacing the bare parameters by their one-loop asymptotic scaling forms 
(in the present case, only the viscosity v and magnetic resistivity vp receive nontrivial corrections at one-loop order, 
so these are the only MHD parameters that get improved to this order). The free correlation function appearing above 



is defined and calculated in Eqs. (A29) and (A27), where further details can be f ound. The frequency integration 
is straightforwardly performed by means of the residue theorem. From Eqs. ( A29] ) and (A27), we immediately see 



that there are two simple poles in the lower complex frequency plane and two simple poles in the upper half plane. 
Closing the contour in the lower half plane and taking the radius of the semicircular contour to infinity yields the 
desired integral Eq. (|7^), and we obtain (at zero-loop order, that is, for the free non- interacting theory) 
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(81) 



The first line holds for arbitrary spatially correlated noise Do(k). We now renormalization group improve this expres- 
sion by replacing the viscosity and magnetic resistivity by their (one-loop) asymptotic scale dependent forms u(k) 
and vsik) calculated above in Eqs. (p4[) and (|55|) to obtain the one-loop spectral density 



£(*) = (— r- ) 



4 



6A d 



(2?r) d V4 + y 

1 1 
1 1 



(A-B)i 



fc-§+|(d-») 



1 -1 
-1 1 



(82) 



This has been evaluated using the specific noise spectrum as written in Eq. ( A29 ). The dimension dependent constant 
Ad is given in Eq. @. From this expression we can immediately read off the individual spectral energy functions 
for each Elsasser variable since each entry of this array Eq. (^|]) involves a quadratic pairing of the Elsasser fields, 
that is 



E(k) 



( Epp(k) 
\E QP (k) 



E PQ (k) 
E QQ (k) 



(83) 



By means of the transformations given in Eqs. ( pl| ) and ( p2| ) we readily obtain the spectral energy densities for the 
physical fluid velocity and magnetic fields, 



E v (k) = {E PP (k) + E PQ (k)} , 
E B (k) = {E PP (k) - E PQ (k)} , 



(84) 
(85) 



respectively. With these expressions in hand, we can now obtain in detail the relative influence of kinetic and magnetic 
noise. From Eqs. (p2) , ( 84 ) , and (Bq), we deduce the individual spectral energy densities, namely 



Ev(k) = ( 

for the kinetic energy spectrum and 



d-l, S d 



(2tt) c 



j+t(d-v)| 



E B (k) - (— - ) 



u s d 



(2ir) d 



(86) 



(87) 



for the magnetic energy density spectrum. The ratio of the magnetic to kinetic energy spectra is given by 

E B {k) fD t 



E v (k) 



demonstrating that this ratio depends only on the ratio of the magnitudes of magnetic and kinetic noises. 

As promised, we have succeeded in calculating the scaling functions for both the fluid velocity and magnetic field 
spectral functions. To one- loop, these results may be read off from Eqs. ( j36| ) and (|87|). Moreover, the difference in 
the two spectral functions is due solely from the dependence on the noise amplitudes. Thus for example, in the limit 
of vanishing magnetic noise, we have A = B. In this limit, the magnetic resistivity Eq. ( |55|) does not renormalize 
and moreover, the magnetic spectral energy density Eq. (|87| ) vanishes identically. Conversely, in the limit of zero 
kinetic noise, the fluid viscosity Eq. ( |54"|) would not renormalize and the kinetic energy density Eq. (|S6| ) would vanish 
identically. Another observation is that when both kinetic and magnetic noises are present, both spectral functions 
scale in the same way. 

The scaling behavior of the energy spectrum has a special interest due to the significance of the Kolmogorov 
spectrum as well as others, such as the Iroshnikov-Kraichnan spectrum. In hydrodynamics, for a freely decaying 
turbulent fluid, which in particular means without an external driving force, Kolmogorov gave general arguments 
that the energy spectrum would behave as fc -5 / 3 . For forced hydrodynamics and/or MHD, it is clear from Eq. ( |8^ ) 
that the energy spectrum can have a range of behavior, with the 5/3 form not being unique. In particular, for a 
specific choice, namely y = d, the Kolmogorov spectrum emerges. However, since the validity of the present RG 
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analysis only is in the asymptotic region, this spectral behavior only occurs in the k — > limit. This differs from the 
region relevant for freely decaying turbulence, which extends from some intermediate value of k up to the higher k 
regime. Nevertheless, an interesting hypothesis can be motivated based on Kolmogorov's universality arguments for 
this spectrum. If the energy spectrum in the forced case has the same behavior as found for freely decaying turbulence, 
then the asymptotic limit of other hydrodynamic quantities, such as the Prandtl number, Reynolds number, and the 
viscosity also should be the same in both cases. Following this hypothesis, it would allow the stochastic hydrodynamic 
equations for y = d to be used as a tool for studying various properties about freely decaying turbulence. Similarly, 
other types of noise (e.g., Iroshnikov-Kraichnan) could be studied, by choosing the noise exponent y accordingly. In 
this paper, this direction of reasoning will not be developed or applied, but it is worth noting it here. 



VI. MAGNETIC HELICITY 



Up to now we have been considering the effects of both random forces and random currents in turbulent MHD 
and the observable consequences these lead to at large spatial separations and for long times. One can also assess 
the influence of magnetic helicity on the asymptotic behavior of the magnetized turbulent fluid. This is especially 
interesting since helicity can arise naturally in rotating turbulent systems such as are encountered in astrophysics 
(e.g., in accretion disks, see Papaloizou and Lin 1995) and galactic dynamics. In MHD, the magnetic helicity is an 
independent quadratic invariant conserved in the limit of nondissipative turbulence. In three dimensions (d = 3), it 
is given by 

H = - [ d 3 x A ■ B, (89) 



2 . 

where A is the vector potential. In 2d MHD, the corresponding quadratic invariant is given instead by the square of 
the magnetic potential (Biskamp 1993). For /i?/drodynamics, the analogous quantity of interest is the kinetic helicity, 
but in nondissipative MHD, this is not a conserved quantity. 

In much the same way as was done for the energy spectral functions derived above, one can derive a spectral density 
for magnetic helicity and use this to study the helicity spectrum at very large spatial and temporal scales. We will do 
so for d = 3 dimensions. First, we note that Eq. ( p9| ) is a gauge-invariant quantity. Thus, by means of the relation 
between magnetic field and vector potential B = V x A, in the Coulomb gauge we can write 

1 /^VxB(y) 



A(x) = — / d A y . (90) 

4tt J \x-y\ 



With this, the helicity can be expressed entirely in terms of the magnetic field as 

1 f^fVxBi 



H = — d x / d 6 y^-^ — ^ • B(x). (91) 
8tt J J y \x-y\ K ' V ; 

By Fourier transforming the magnetic field and the vector potential Eq. (^) to wavenumber space, it allows the 
helicity to be expressed as 

1 f d 3 k -. _ -> -. 

H = ^ J j^k- 2 [ikxB(k,t)}.B(-k,t). (92) 

We now define a helicity density (helicity per unit spatial volume) and an associated helicity spectral density. For 
this, first Fourier transform the integrand in Eq. (|92|) to frequency variables and take the stochastic average of H 
over the random fluctuations. This allows the averaged magnetic helicity to be expressed as 

/>oo 

H = V dk H(k), (93) 
Jo 

where V = (27r) 3 <5 3 (0) is once again a spatial volume factor and the helicity spectral density is given by 

r ./( ) poo 7 

H(k) = 8^ J (2^ loo 2^ e ^ MCp ™ p > ^ W) " Gp -^ ^ (94) 

Here, £y m is the fully antisymmetric tensor density in three dimensions. The indicated integrations are taken over 
the unit sphere and frequency, respectively. In arriving at this final expression, we have averaged over all the sources 
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of random fluctuations. We also point out that H is a real quantity; this can be checked directly in Eq. (|||) by 
using B(k, t)* = B(—k, t) together with the symmetry properties of tij m - In arriving at this final expression, we have 
written the magnetic field correlation function in terms of Elsasser field correlations (Cpp, Cpq). 

It is clear that the magnetic helicity spectral density (as well as the total helicity) vanishes identically in a system 
subject strictly to non-helical forces and currents. This is because the resulting magnetic field correlation function 
is proportional to a symmetric tensor, < B m (k,t)Bi(—k,t) > oc P rn i(k). The same symmetry characteristics are of 
course shared by the correlation functions for the Elsasser variables. It is furthermore clear from Eq. (|94|) that to have 
a net averaged helicity, the correlation function must contain an antisymmetric contribution, and that this helical 
contribution must in fact be proportional to the epsilon tensor. In other words, to include helicity one must "build" 
an antisymmetric contribution to the correlation function using only the the metric, or Kronecker delta <5y, products 

of the wavevector k m , and epsilon itself ei n k- This we will do dynamically by injecting an initial spectrum of random 
magnetic helicity into the MHD equations by means of a straightforward extension of the noise term. We will then 
track this random helicity dynamically by means of the renormalization group applied to helical MHD and compute 
the asymptotic limit of the helicity spectrum. 

The presenc e of random helicity in three dimensions can be incorporated in the system by extending the noise 
spectrum Eq. ( A25| ) as 



[T (k, uj)] mn = D Q {k) Pmn(k) - iF (k)e mn i ki, (95) 



where the nonhelical contribution arising from the random forces and currents is that given in Eq. ( |A25| ) and the 



second term involving the epsilon tensor is due entirely to helicity. In higher (d > 3) spatial dimensions, the "extra" 
indices on the corresponding epsilon tensor would b e con tracted by additional factors of the momentum vector. The 



noise function (matrix) -Do(fc) is specified in Eq. (A29). The noise function (matrix) associated with the helical 



contribution can be modeled by taking the following initial power law spectrum for random helicity, 



F (k) =2h B k~ 



1 -1 
-1 1 



(96) 



where hp > is the amplitude of the helical noise and the number w is an exponent (analogous to the exponent y 
characterizing nonhelical noise). Both are free parameters. We have written the helical noise directly in the Elsasser 
field basis Eq. (|TT|). In terms of the physical fields, this source Eq. (^) represents adding a helical noise term to 
the right hand side of Eq. (0), leaving both Eqs. (||) and (fjj) unaltered. One can add a random kinetic helicity to 
the right hand side of Eq. (j^) if one wishes, by a simple modification of Eq. (|96|). To our knowledge, in contrast to 
the kinds of scaling arguments leading to the derivation of the Kolmogorov spectrum, there are no arguments that 
can shed light on what kind of helicity spectrum one might expect to see for example in a freely decaying turbulent 
helical fluid. As such, there are no specifically interesting values of w of which to take note. 

In a series of Appendices, we have set up and carried out the renormalization of the response, noise and vertex 
functions associated with the MHD equations in the absence of helicity, Kb = 0. Now that we have modified the noise 
spectrum in Eq. (^5|), we must repeat the renormalization group calculations taking into account the new helical 
contribution. If we refer to the elementary Feynman diagrams in Fig. (la)-(lc), we immediately see that the only 
elementary diagram to be modified is the correlation function, Fig (lb). On the other hand, the bare response and 
verte x functions are unaffected by the presence of Eq. (|96|). When computing the loop diagrams, rather than Eq. 
( A29| ), one must use 



[Co (k,u>)] mn = C {k,u))P mn (k) - iCg(k,uj)e mn ik h (97) 



where Cq is obtained from Cq in Eq. (A29) by replacing Do(k) by Fo(k). In terms of the diagrammatic language, 
the response function diagram, Fig. (2a), and vertex function diagram, Fig. (2c), are unaffected, while the "blob" 
in the diagram Fig. (2b) for the correlation function is modified by the addition (insertion) of the helical noise. All 
the calculations in the Appendices are to be carried out with this modified diagram and setting d — 3, as we are 
considering helical noise in three dimensions. The details and mechanics of the calculations are very similar to those 
presented in the Appendices. The net result of this procedure simply is that at one loop order the random helicity in 
Eq. ( p5| ) does not renormalize the response function, the noise spectrum nor the vertex function in the hydrodynamic 
limit, i.e., for k — > and for t — > oo. In the case of the response and vertex functions, the angular integrations are 
responsible for yielding null results in the hydrodynamic limit. The non- helical noise spectrum does not renormalize 
for w > — 1, which is a mild condition we impose on this exponent. This means that the original set of renormalization 
group equations Eq. ( [f0| ) (as well as their solutions, fixed points and stability properties) remain unaffected (at least 
at one-loop) by the presence of random magnetic helicity. 
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Independently from this, observe that the helicity density never contributes to the system's total energy Eq. ( J7q ) 
(since helicity is traceless). Conversely, the energy density never contributes to the helicity (the energy is linked to a 
purely symmetric tensor). Nevertheless, the helicity has an impact on the dynamics of the system, generating its own 
contribution to the velocity and magnetic field correlations, as is clear from Eq. (|97|). We come back to this important 
point after calculating the helicity spectral density by means of RG-improved perturbation theory at one-loop. 



A. Helicity spectral density 



In analogy with the energy spectrum, we can calculate the helicity spectrum for both, non-interacting and interacting 
dynamics. To this end, we define the following helicity (matrix) function (in three dimensions) 



V 



d 3 k 

W) 

= V dkn{k) 
Jo 



dto 



3 ' 2-7T 



. k 3 [C (k, 



where the helicity spectral density function (also a 2 x 2 matrix) is given by 



H(k) = 



dtt 3 f°° du 



(98) 



(99) 



We need to know the correlation function (for the Elsasser fields) [C (k,ui)] m j in order to determine the helicity 
spectrum. We will obtain this spectrum by means of (1-loop) RG-improved perturbation theory. So, we first evaluate 
the spectrum Eq. ( p9| ) in the non-interacting field limit and then improve the resultant expression in order to obtain 

the corrected spectrum, which is valid at one-loop. For the first step, we simply replace [C {k,u))] m j by its free-field 
limit as given in Eq. (|97|). Note in comparing Eq. (97) with Eqs. ( |98|) and (|99|), we see that only the antisymmetric 
part of the full correlation function contributes to the helicity and to the helicity spectrum, which should come as no 
surprise. Carrying out this replacement yields for the zeroth-order helicity spectrum 



Ho (k) 



s. 



duJ rFtt \ 
Z7T 



47r(27r) 3 

where the factor 63 results from the integration over angles (in d = 3), and we must integrate 

C F (fc,w) = g (k 1 u;)F (k)g (-k, -u) 



(100) 



(101) 



over the frequency. This involves the identical frequency integration performed for the calculation of the energy 
spectrum, in arriving at Eq. (81). Thus, direct use can be made of that step, without any further calculation, to 
obtain 
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(102) 



for the zeroth-order value of the spectrum (which is also the spectrum corresponding to the trivial fixed point) . The 
matrix products are worked out making use of Eq. (Bq) . The RG-improved helicity spectrum is therefore given by 



n PQ (k)\ 

h B S 3 (6A 3 (A-B) 



4tt (2tt) 



— ( 

7r) 3 V 



1 



y 



1 1 



(103) 



which indicates this helicity spectral matrix refers to the Elsasser P, Q field basis. In terms of the physical fields, the 
magnetic helicity spectral density is given by 
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H M (k) =Hp P (k)-H PQ (k) 

hB_S3_/ 6A 3 D B yi + i {1+y) 

2ir (2^) 3 V 1 + y J ' ( ' 

We note that the other linear combination yields the kinetic helicity spectrum, but this is identically zero, Hy — 
Tipp +H-PQ = 0. This is consistent with the fact that no kinetic helicity has been injected into the dynamics. We see 
the difference in scaling of the helical spectrum as determined from the fixed points: apart from numerical prefactors, 
the magnetic helicity scales as 

n M (k) ~ k~ w , (105) 
for the free, non-interacting theory (trivial fixed point) while it scales as 

n M (k) ~k~ w+ ^ 1+ y\ (106) 

at the IR stable nontrivial fixed point (interacting theory). 



B. Renormalization Group Improved Correlation Functions 



Just as the general scaling analysis used for obtaining the general form of the energy density spectra is rigorously 
complemented by improved perturbation theory, the scaling behavior derived above for the field correlation functions 
can be complemented by renormalization group improvement. Thus, the scaling form of the correlations displayed 
in Eqs. fl63| ) and ( |64[ ) is correct, as it stands. However, this only makes explicit the scaling exponent, while lumping 
together with the scaling function the as yet unknown proportionality constants. To obtain the latter, a bit more 
work is required. Once again, renormalization group improved perturbation theory leads directly to the answer. 

Having introduced the magnetic helicity, we opt for improving the full correlation function Eq. , which contains 
both kinetic and helicity effects, since the same basic steps are needed for improving each contribution separately. 
Improvement of a particular quantity consists, as we have already seen, in replacing the bare viscosity and magnetic 
resistivity by their (one- loop) running expressions Eqs. (54) and ([So]). 

Fro m the structure in Eq. (|97|), we see that this is tantamount to improving the produ ct of response functions Eq. 



(|A27|) that appear in both the Cq and Cq contributions, as can be checked in Eq. (A2S). Since both noise matrices 



commute with the product of response function matrices, the problem reduces to improving the following matrix 
function 



g (k, u)Qo(-k, -u) = f(v, v B ; k, u) 

1+k- 4 



1 
1 



f{y,v B ;k,u) = ■ ^ Y2 C \T2h 

(7+ + T-) (7+-7-r Kk J 

where the dimensionless (in the free-field limit) scaling function C is explicitly given by 



(107) 
(108) 
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\k 2 
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[fe 2 (7++7-)l 2 



)(1 



l^(7+-7-)l 1 



(109) 



We make no special distinction between functions of 7+, 7_ and functions of v, vp, since these viscosity and resistivity 
parameters a re lin early related via Eq. (^2|), and are freely interchangeable. 

From Eq. (107) and the commutativity of the block matrices, we can express Eq. (J97|) as 



[Co (*,w)] T 



(PP) (PQ) 

(QP) (QQ) 



= f{v,v B] k,u)\2k- y 



A B 
B A 



Pmn(fc) + 2h B k~ 



1 -1 
-1 1 



(110) 



Of course, the term proportional to hp holds only in d — 3, while the remainder of this expression can be evaluated 
in any spatial dimension. This yields the individual correlations in terms of the Elsasser fields, but we can skip 
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immediately to the quadratic correlations involving the fluid velocity and magnetic fields by means of the linear 
relations given in Eq. (fill). From this, it implies 



(v lVj ) = f{y, v B \ k, to) (A + B)k-y (k) + 



I B 



{BiB-) = f(v, v B \ k, u) (A - B)k-y P„ (k) 



(viBj) = f(v, u B ; k, to) Hb k~ 



A + B 

h B 
A-B 



k y - 
k y - 



(in) 

(112) 
(113) 



The improved correlations are obtained immediately upon making the replacements mentioned above. We see that 
this operation affects only the overall factor f(y,VB]k,oS), which is common to each individu al co rr elati on func tion . 
As such, this is the only function that needs to be improved, with the correlations in Eqs. (Ill), (112), and (115) 
improved by replacing 



From Eqs. (107), (M) and 



we find that 



(114) 



f(v(k),v B (k);k,u) = 
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where the constants 



>6A d (A + B) 
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6A d (A - B) 



(115) 



(116) 



We point out t hat t he scaling function C is also modified by the running of the viscosity and resistivity as should be 
clear from Eq. ( 109 ) , and this is acc ount ed for by the expo nent z. 

Substituting Eq. (|115|) into Eqs. (Ill), (112), and (113) gives the complete structure of the correlation func tions 



corre ct at one-loop order. Note that the net exponent —4 + 4(4 + y — d) — y = —d — 2\ — z, so that Eqs. (Ill) and 



( p.12 ) reproduce the same scaling (for h B = 0, i.e. zero helicity) as we found in Eqs. (]6Sj) and (|64j). However, we 
also explicitly obtain the associated scaling functions. We see that the presence of random magnetic helicity modifies 
these correlations, since the helical contribution scales as k v ~ w+1 relative to the non-helical part. Furthermore, 
the relative strength of this term goes as the ratio of the amplitudes of the helical to the non-helical fluctuations: 
h B /(A + B),h B /{A-B). 



VII. THE CALLAN-SYMANZIK EQUATION 



In this section, we derive the scaling behavior of the correlation functions in the Fourier domain from a distinct but 
entirely complementary point-of-view provided by the methods of field theory. In this approach, the renormalization 
of the ultraviolet (short-distance) divergences of the stochastic field theory for MHD formulated in Appendix A leads 
to certain partial differential equations, frequently referred to as Callan-Symanzik equations^ for the one-particle 
irreducible (1PI) Green functions. The exact solutions of these equations provide considerable information that goes 
well beyond perturbation theory. These solutions tell us how the 1PI functions scale near the fixed points of the 
theory and permit one to quantify the approach to the fixed points. Moreover, the concept and mechanics of RG- 
improvement is given a rigorous basis by means of the Callan-Symanzik equation. We therefore devote some space 



1 Equations of this type are typically referred to as Callan-Symanzik equations, though strictly speaking, this name refers to 
a specific type of renormalization group equation resulting from so-called mass-subtraction (in quantum field theory). A mass- 
independent subtraction leads to a different equation, similar to the one actually derived here, and might be more properly 
called a t'Hooft- Weinberg equation. A clear distinction between Gell-Mann-Low, Callan-Symanzik and t'Hooft- Weinberg 
renormalization group equations is drawn by (Gross 1981). 
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in the present section to an application of these ideas to stochastic MHD field theory. In what follows, we consider 
non-helical MHD. The extension to the helical case is straightforward but will not be carried out here. 

To begin, we first tabulate the canonical dimensions (in units of momentum) of all the parameters and fields 
appearing in the generating functional for MHD. These are easily obtained by noting that the action A in the 
dynamic functional, 

Z = J [Va][V$\ cxp(A), (117) 

must be dimensionless: [A] — 0. We recall that we have employed a short-hand notation for the pair of Elsasser fields: 
<J? = (P, Q), and the conjugate fields a arc defined in Appendix A. We can read off the explicit structure of the action 



A directly from Eqs. (A14), (A18), (A19) and (A25). From this zero-dimensionality condition we therefore find that 



(where [X] = dx denotes the dimension of X in units of momentum) 

h±]=cLy, [w]=tL = 2 + d 7 , (118) 

[A} = [B}=d A = d B , [T]=d r = d A -y, (119) 

[A] = d x = (£z|±l!) + - \d A , [g ± }=0. (120) 
The Elsasser and conjugate fields have dimensions 



[*] = d» = -3 - | - ^ + id r , [8]=d„ = -l-i- ^ - Jdr, ( 121 i 



respectively. Note that the canonical dimension of the fields depends on the noise exponent y since dr = d A — V- 
Finally, the two 1PI functions of particular interest to us are the response function and the noise spectral function. 
The response function has dimension 

Fn] = rfn = [(Go \j] = (2 + Ay), (122) 
while the non-helical part of the noise spectral function has dimension 

Fra] = d 2 = [{T)ij} =d A -y. (123) 

We remark that it is customary to denote the 1PI Green functions by a capital T, which is not to be confused with 
the block array of noise correlation functions. These 1PI functions are built up out of ce rtain numbers of conjugate 



fields and physical fields: the indices on the first two terms on l eft hand side of Eq. (122) indicate that one (1) 



conjugate and one (1) physical field are involved, while in Eq. (123), no conjugate field (0) and two (2) physical fields 
are involved. We suppress the latin indices here to avoid clutter. Note all canonical dimensions can be expressed in 
terms of the dimension of the viscosity (and resistivity) d 7 , the noise amplitude d A = ds, and exponent y. 

The one-loop corrections to the MHD parameters calculated in the Appendices have been cut-off with a momentum 
or wavenumber regulator A. Recall, this actually is a physical cut-off in momentum. It is associated with the 
smallest length scale below which the continuum description of the fluid breaks down and must be replaced by an 
atomic/molecular description. If this short-distance cutoff is taken to zero (A -1 — * 0), i.e. take A — * oo, it would 
result in short distance or ultraviolet divergences showing up in the one-loop response and correlation functions that 
are calculated in the Appendices. To systematically remove these ultraviolet divergences, renormalization constants 
would have to be introduced in the parameters and fields that appear in the above action A. This may be implemented 
by introducing renormalized parameters (denoted by Z a , Z$, Z a ) and renormalized fields (denoted by an R) as (Gross 
1981; Amit 1978) 

ta = Z- x £, (124) 
* = 4 /2 * fl , (125) 
a = Zl' 2 a R . (126) 

Here the symbol £ a is a shorthand notation denoting the collection of all parameters that appear in the dynamical 
equations, £ a = (A, P,7_|_,7_, A), with the label running from a — 1,..,5. We will be able to draw some general 
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conclusions about the asymptotic properties of the correlation functions without having to explicitly calculate these 
renormalization parameters. 

In the following, we first consider an arbitrary 1PI function, and then specialize to our response and noise functions 
to obtain information regarding the correlation function. Keep in mind that for the MHD theory considered here, 
these functions all have the structure of two-by-two block arrays multiplied by a d-by-d tensor factor. We suppress 
the tensor indices in what follows. The renormalization of the 1PI functions containing N factors of the physical field 
and N factors of the conjugate field is given by 

T B n ,n (k, w; g; n) = Z N J 2 Z§' 2 T N ^ (k, w; {&}, So; A). (127) 
In arriving at this expression, we have expressed the bare fields in terms of the renormalized fields by means of 



Eqs. ( Jl 25| ) and (126). The label R denotes a renormalized quantity. Most importantly, the renormalization requires 
introducing a finite but arbitrary momentum scale fi, which remains in the finite, renormalized expressions. As it 
stands, this equation merely expresses the fact that the function in question is renormalizable. In other words, the 
UV divergences can be consistently subtracted out. This subtraction is done at some arbitrarily chosen scale fj,, and 
this fact must be recorded in the resulting finite renormalized function. Using the fact that the bare, unrenormalized 
1PI function does not depend on this arbitrary momentum scale /x (it depends only on the ultraviolet cut-off A), we 
can easily derive a differential equation for the renormalized 1PI function as follows: 



>^) o T N ^ {k,u;{£ a },g ;A) = 0. (128) 



This implies (see footnote 1) 



+ P{g) Tg + ^ Aa ^ 9^ " T 7 " f 7 ) ^ ^ W; { ^ } ' 9; M) = °' (129) 
where the coefficient functions appearing in this differential equation are defined by 

A« = (M^) . (130) 

m = (m^) . (131) 

, <91nZ$, „ , d In Z a > 

7= ( M ^H°' ^=^-djr)o- ( 132 ) 

The zero-subscript (0) indicates the derivatives are taken while holding fixed the bare parameters This equation Eq. 



(129) indicates that a change in the arbitra ry b ut finite subtraction scale [i is compensated for by corresponding 
changes in the physical parameters. As Eq. ( |l29| ) is a linear first-order partial differential equation, it can be readily 
solved by means of characteristics (Amit 1978). The exact solution is given by 



^(M; M) = e*p(-y — [y7(«) + y7W]J T n ,n (k,uj;{^{b)},g(b);bn), (133) 
where the equations for the characteristics are 

P(9(b)) = b^, (134) 
dln£ R (b) 

A a (g(b)) = b am ^ 0) . (135) 

These equations are subject to the bou ndar y conditions at b = 1, g(l) — g and ^(1) = £,a ■ The solution can be 



readily verified by differentiating Eq. (133) on both sides with respect to b and using the chain rule. The set of 

characteristics define a ruled (hyper) surface on which r N ^ is defined and b parametrizes the characteristic curves 
on this surface. It is related to the RG flow parameter via b = e~ e . This dimensionless parameter b can be related 



to ratios of momentum (or length) scales, as we will see below. The physical content of the solution Eq. (133) is 
as follows. A change in momentum scale (as measured by b) is accompanied by a change in the coupling g(b) and 
in the other parameters (b). Furthermore, the vertex function picks up an additional overall exponential factor 
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depending on the two anomalous dimensions 7 and 7 that are defined and calculated in Eq. (132). In particular, we 
see that the anomalous dimensions depend on the renormalization of the physical and conjugate fields (wavefunction 
renormalization) . These dimensions control in turn, the approach to the fixed point, as we will see below. In a nutshell, 
this equation relates the 1PI function at two different scales. We now use this fact to deduce how the correlation 
functions depend on scale. 

Define dimensionless scaling functions {d N pj is the canonical dimension of the 1PI function) F NN as 



<-»-R 



T NN {k,w;{g},g;ii)= fi° 



F 



NN 



2+d-, 



/ 



(136) 



so that the dimension of IV at * s carr i e d entirely by an appropriate power of /z. The individual arguments of F N fi 
are rendered d imen sionless since each is expressed as a dimensionless quotient, by making use of the dimensions listed 
above in Eq. (118). Here, d^ a denotes one of cIa, rfs, d\, d 1+ , <£y_ • Then from the solution Eq. (133) 



T N . N (fc,w;{£f}, ff ;^) = (bfx) 



exp 1 



x Fnn I 7— , 



dx.N . , N_, 

— hr7W + tt7W 
! x 2 2 

r ffW 1 m 



(137) 



At this point, we have succeeded in writing each renormalized 1PI Green function as a product involving a dimcn- 
sionful power of the momentum scale fi times an exponential factor times a dimensionless scaling function. Now, the 
correlation function for the Elsasser fields P, Q, which recall is not 1PI, is given by (see Eq. (A29)) the block-matrix 

product C= (Tn)- 1 F02 (f Thus, from Eq. pg) 



C(k,iu;{^},g; f i) = (b f i) 



dA-y-Zd-, 



4 «p( + ^§ 7(*))<?(^, 



2+d. 



e(6) 

(6/i)<^ 



};s(6) 



(138) 



where C = (Fn) _1 F02 (-Fn) -1 and from Eqs. (122) and (123) d 02 - 2dn = d A - y - 4 - 2d 1 . Recall that C is a 



2x2 block array with a pair of <i-dimensional indices. Note that in Eq. (138) the dependence on the conjugate field 
anomalous dimension 7 has dropped out. 



It is worth commenting that provided the characteristic equations Eqs. (134) and (135) are solved, Eq. (138) yields 
the exact solution for the correlation functions in terms of the parameter b which is a measure of the momentum 
scale (or inverse length scale) at which one observes the system. The exact solution is known if explicit knowledge is 
available of the anomalous dimension j(x) and the solutions of the above characteristic equations. However, typically 
these are known only up to some order in perturbation theory, which usually is to a given number of loops. Thus, 
the solution of the Callan-Symanzik equation is only known to the same number of loops. The important utility of 
the solution lies in the fact that it tells us how quantities behave away from the fixed point and it places the results 
obtained from naive scaling arguments on a solid footing. In general, the beta function Eq. (131) will not vanis h fo r 
arbitrary b (i.e., far away from fixed points) nor will the coupling g(b) be constant. As such, the expression Eq. (138) 
tells how the fields are correlated in general, that is, everywhere in parameter space, and not just in the vicinity of 
the fixed points. Of course, if the system happens to be at a fixed point, g — * g* = const., 7^7*= const., etc. 
Choosing b = k/fj, and restoring the block labels and tensor indices gives for the PP entry of the block array 



C 



PiPr. 



(k 1 u J ) = k d --y- 2d ^ 4 (-) p jn (fc)c>p( 



UJ 



and similarly for the PQ and QQ entries. At a fixed point, i.e. A* = A a (g*), 



k N A* 



£?(*0 (7) 



(139) 



(140) 



which follows from evaluating Eq. (135) at a fixed point. In particular we see that the correlation function, calculated 
up to n-loops, is given by the free correlation function but expressed in terms of the scale dependent (running) 
couplings and parameters (calculated up to n-loop order). This is nothing but the prescription leading to RG-improved 
quantities, but her e it is a natural and automatic property of the solutions of the Callan-Symanzik equation. 
Comparing (13£) to Eqs. ( |63|) and ( |64| ) yields the exponent relations 
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z = 2 + d, 



11 



i d l t, 



d,A + y--f*)- 



(141) 



These relations allow the various "engineering" dimensions to be related with the fixed point exponents calculated 
earlier. Moreover, we see how the concept of RG-improvement used throug hou t this paper appears as a rigorous feature 
of the solutions of the Callan-Symanzik equation. The equality in Eq. (139) shows that the correlation function in 



the vicinity of the fixed point is a power law multiplied by the scaling function. Furthermore, the scaling function is 
evaluated with the running or scale-depende nt va lues of the parameters, that are computed in the neighborhood of 
that fixed point, as exactly indicated by Eq. ( |140| ). Note that if there is any nontrivial wavefun ction renormalization 
(indicated by a 7* ^ 0), it is automatically subsumed in the exponent \ as shown in Eq. (141). Since we have 
calculated both z and x Dv independent means, there is no need for a separate calculation of 7* if we are interested 
only in what is happening at the fixed point. 



A. Approach to the fixed point 

Observe that the exponent (also called the anomalous dimension) j(g) controls the approach to the scaling regime, 
which is reflected in the pure power law behavior of the correlation function at the fixed point. In general, it is 
unlikely the dynamical system will be sitting exactly at the fixed point. Therefore one should study the form of the 
correlation function for values of the coupling that arc in the basin of attraction of the fixed point and see how this 
function depends on arbitrary initial values of the coupling. 



This can be done systematically by returning to the exponential factor in Eq. (138) and carefully expanding this 
about the (non-trivial) fixed point. First make a convenient change of variables as 

-j(g(x))= / dl 7 (g(l))= / jgldy, (142) 

X Jo Jg(0) P{y) 

where b = e l . Note, this I is not the same as the RG flow parameter I in Eq. (ftcf) , but they are related as — ln(6) = I 
(Frey and Tauber 1994). As already pointed out, b parametrizes the sense of "flow" along the characteristic curves 
that arise in the solution of the CS equation and does not automatically single out either the infrared (decreasing b) 
or the ultraviolet (increasing b) limit. The connection to these physical limits is made once we assign infrared and 
ultraviolet flow "directions" to the characteristics. The one-loop beta function can be read off directly from Eq. ( |4^ ) 
and is given by (3(g) = I2g(g* — g). Solving for the coupling by integrating the differential equation Eq. ( f43| ) gives 

as. . ^ exp(rf) . (Ms, 

where e = 4 — d + y > 0. The infrared stable non-trivial fixed point is reached by taking b = ^ — > or I — > —00. In 
terms of the Wilsonian flow parameter £, this corresponds to the limit I — > 00. Note this fixed point is reached for all 
initial values g(0) of the coupling. 

Next, Taylor expand the anomalous dimension about the fixed point as 

7(5) = 7* + 7' (9 ~ 9*) + 0(g - g*) 2 - (144) 



With this, evaluating the integral 



■"TW^-tUm, ,145, 



>,m m " 1 12 y S (o) 

it implies the exponential factor contributing to the correlation function Eq. ( |l38| ) works out to be 



exp 



7 



(146) 



Recall that the fc-dependent factor with ex pone nt 7* is already subsumed into the definition of the fixed point 
"roughness" exponent x, as indicated in Eq. fll4l|). However the second term, which depends on the exponents 7' and 
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e, clearly is a correction to the fixed point (power law) form of the correlation function written in Eq. (139). Depending 
on the algebraic sign of 7', this extra term, which controls the approach to the fixed point, can either augment or 
diminish the energy spectrum with respect to its asymptotic fixed point value Eq. (|7q). Similar consideration holds 
at the nontrivial fixed point for the asymptotic form of the helicity spectrum in Eq. ( |106[) . 

Some comments are in order about the region away from the fixed point. First, any amplification or diminution 
in the spectra densities with respect to their steady state fixed point values, which is brought about by this extra 
term, is a statement about how the densities evolve as they approach the nontrivial fixed point. It is a conceivable 
possibility that numerical simulations, due to their spatial and temporal limitations, actually may be seeing some 
sufficiently long-lived non-steady state of this type 



This point will be elaborated upon in Sec. IX 



However, we 

hasten to point out that if the wavefunction renormalization constant turns out to be equal to unity, Z$ = 1 , then 
t he a nomalous dimension Eq. (132) is identically zero, 7 = 0. When this is the case, then the expontial factor in 
(138) reduces to unity, and any corrections to pure power law scaling will be due entirely to the scaling function itself. 
For incompressible hydrodynamics subject to random forces, it is known that Z = 1 (DcDominicis and Martin 1979). 
This in turn can be shown to be a consequence of the underlying Galilean invariance of the Navier-Stokes equation, 
which is maintained for all spatially correlated random forces. Si nce t he equations of nonrelativistic MHD are also 
Galilean invariant, we suspect that the anomalous dimension Eq. (132) will also be vanishing. This issue will not be 
discussed further in this paper. 



VIII. CASCADE DIRECTIONS 



One important question of general interest in the study of hydrodynamic systems is how the nonlinear terms 
redistribute energy and, if present, helicity amongst the various scales. In particular, if energy and/or helicity is 
introduced into the system at a given scale, it then will be distributed by the subsequent evolution of the hydrodynamic 
or MHD equations, either to smaller (direct or normal cascade) or larger (inverse cascade) length scales. In general, 
this discussion of cascades is restricted to the inertial range, which corresponds to length scales sufficiently large that 
they are not significantly dissipated by the viscosity term. 

It is interesting to reflect on the history of the subject of cascades. The earliest expectation, based on simulations 
and general arguments of the Navier-Stokes equation in d = 3, was that a direct cascade (also suggestively called 
a normal cascade) would occur from hydrodynamic evolution. Subsequently in (Kraichnan 1967), for the Navier- 
Stokes equation in d — 2, based on general observations about the two conserved constants of motion, cnstrophy 
and kinetic energy per unit mass, it was argued that an inverse cascade of energy should occur. For MHD, these 
arguments of inverse cascade are not applicable in either d = 2,3 and so a direct cascade of energy is expected 
(Biskamp 1993; Kraichnan 1973). On the other hand for MHD, considerations based on conservation of magnetic 
helicity in d = 3 and its counterpart, the square of the magnetic potential in d = 2 showed that an inverse cascade 
is expected of the respective quantities (Frisch et al. 1975). All the above expectations are based on compelling, 
but not rigorous, arguments. Verification of these expectations mainly have come from computer simulations directly 
of the hydrodynamic/MHD equations or indirectly through some approximations to these equations (Brandenburg 
et al. 1996; Lilly 1969; Frisch and Sulem 1984; Sommeria 1986; Pouquet et al. 1976; Pouquet and Patterson 1978; 
Pouquet 1978). Clearly this approach has its limitations in the interesting regime of high Reynolds number, due to 
finite simulation time and system size restrictions. 

Returning from this digression back to our RG analysis, cascades clearly are associated with time dependent 
processes. On the other hand, the RG analysis in this paper examined steady state solutions of MHD with time 
translationally invariant force terms and only treated the leading asymptotic behavior at large distance. From this 
analysis, information about the energy spectrum can be obtained for the largest length scale. In particular, this means 
length scales beyond the largest characteristic scale at which energy is input from the external force. As such, the 
RG analysis of this paper can determine presence or absence of inverse cascade, depending whether the infra-red limit 
of the energy spectrum is enhanced or not. This is a useful application for RG, since generally some direct cascade 
always occurs and the interesting question is whether there also is some inverse cascade. 

In order to extract information about cascades from this framework, a few points must be appreciated. The force 
term is continuosly inputing energy and/or helicity into the system with some prescribed spectrum. A cascade for 
such a case therefore is associated with the difference in the steady state magnetic and velocity spectra directly created 
with the force (in the absence of mode mixing) versus the steady state spectra that emerges from the subsequent 
MHD evolution. For example, if the spectra associated with the force is less enhanced in the ultraviolet versus the 
resulting spectra after MHD evolution, then the MHD dynamics is directing energy flow to smaller length scales thus 
exhibiting a direct cascade. Similarily, if the input spectra associated with the force is less enhanced in the infrared 
versus the resulting spectra from MHD, then MHD is exhibiting an inverse cascade. 
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To proceed with this approach, the implied magnetic and velocity spectra associated with the force terms must 
be determined. The point here is Eqs. (g) - (0) give the spectra of the force, but what we seek is the associated 
spectra this implies for the magnetic and velocity fields. To determine these, MHD with the force terms, but without 
the nonlinear terms, must be examined. These linearized MHD equations express the corresponding magnetic and 
velocity spectra for given input force spectrum, before there is any mode mixing. Equivalently stated, these are the 
associated input magnetic and velocity spectra for the given force, prior to any nonlinear MHD evolution. 

The procedure then for determining information about cascades is to compare for a given force term the corre- 
sponding magnetic and velocity spectra from the linearized MHD equations versus those from the full nonlinear MHD 
equations. Depending on how the ultraviolet and infrared behavior of the latter changes compared to the former, 
cascade directions can be deduced. 

Following this criteria, let us determine the cascade directions predicted at the largest length scales from our RG 
analysis. First note that the spectra for linearized MHD equivalently stated is the one also for the trivial fixed point. 
Therefore one fact immediately follows from this criteria. In the region of stability of the trivial fixed point, the RG 
analysis predicts no cascade whatsoever, at the largest length scales. 

Turning next to energy spectra in MHD, from Eqs. (^) and (f75| ) E(k) ~ k a with aiinear = —3 + d — y and 
anontriviai = -5/3 + (2/3) (d - y), so that aiin ear < (>)a n ontriviai in the region 4 - d + y > (<)0. Thus we see that in 
the stability region of the nontrivial fixed point, the input energy spectra ~ P llnc " is more enhanced in the infrared 
and less in the ultraviolet versus the resulting spectra from MHD evolution fc Qnontlivial . Thus our RG analysis of MHD 
predicts direct energy cascades at the largest length scales. Since the above case involves the entire stability region of 
the nontrivial fixed point, this analysis finds at the largest length scales, no prediction of an inverse cascade of energy 
in MHD for any d. 

The RG calculations in this paper easily can be reduced from MHD to hydrodynamics by identifying the fields 
P = Q. The outcome of this is the predictions of the fixed points for the hydrodynamic case are found to be exactly 
the same as given in Sec. IV. As such, our analysis also expects direct but no inverse energy cascades at the largest 
length scales in hydrodynamics for any d. 

It is interesting to note that the change from direct cascade to inverse cascade in the RG analysis occurs right at 
the boundary of stability for the nontrivial fixed point. In particular the unstable region of the nontrivial fixed point 
would predict inverse cascade. In the Conclusion we will discuss interpretations of this observa tion. 

The same analysis as above also can be done for helicity. The helicity spectra give n in Subsect. VI A are TLm (k) ~ k@ 
with /3ii ncar = — w from Eq. (105) and Aiontrivial = — w + (1/3) (1 + y) from Eq. (106). The region of direct helicity 
cascade is /3ii ne ar < /^nontrivial, which implies 1 + y > 0. Since this is for d = 3, note this last condition also corresponds 
to the stability region of the nontrivial fixed point, i.e., 4— d+y = So once again, RG predicts a direct magnetic 

helicity cascade in the entire stability region of the nontrivial fixed point and no cascade in the stability region of the 
trivial fixed point. Furthermore, once again the change from absence to presence of inverse cascade for the nontrivial 
fixed point is at its boundary of stability. 



IX. CONCLUSION 



In this paper a Renormalization Group (RG) analysis has been performed for MHD with a stochastic force term. 
Such an analysis enables the study of highly turbulent MHD in the asymptotic regime of large length and time 
scales. In the format of the RG approach, the central problem is to determine the effective MHD equations at large 
scales that emerge once all short distance physics is integrated out of the initial fundamental MHD equations. The 
analysis was carried out for a general class of force functions Eq. (0) - (Q) and for arbitrary spatial dimension d. The 
outcome of this analysis, as detailed in Sec. IV, showed that MHD for this class of force functions, has two types of 

In one regime, the system becomes a free theory 

and in the other 
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asymptotic behavior, depending on d and the nature of the force. 

in the asymptotic or hydrodynamic limit which is represented by a trivial fixed point in Sec. 
regime the system is a nonlinear interacting theory at the largest scales, w hich is represented by the nontrivial fixed 
point in Sec. IV. For both regimes, correlat ions were computed in Subsec. VB of the magnetic and velocity fields at 



large length and time scales and in Subsec. V C of the en ergy spectrum . In addition, the viscosity coefficients in the 
effective long distance theory were computed in Subsec. VA. Helical force also was treated in Sec . N% Up to the 
one-loop analysis in this paper, the fixed point struc ture w as unchanged with respect to that in Sec. [IVjlThc helicity 
spectrum at both fixed points was given in Subsec. VIA. It would be interesting to test the RG predictions against 
simulations for all these quantities. 

Another interesting direction examin ed w as mixing between velocity and magnetic energy at large scales. The 
results of the one-loop analysis in Sec. V C indicate that in fact at the largest scales magnetic and velocity energy 



will not mix. Thus, for example, if only the velocity force was nonzero, then no magnetic energy would be created at 
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the largest scales and visa- versa. These results of course are only concerned with the asymptotic large scales. Clearly 
the nonlinear interactions in the MHD equations generally should mix velocity and magnetic energy and this is seen 
in simulations (Pouquet et al. 1976; Pouquet and Patterson 1978), which suggest an approximate equipartition of 
energy. However, the RG analysis indicates that at increasingly larger length scales this mixing diminishes. It would 
be interesting to see if such trends also are seen in simulations. 

The RG analysis of this paper address ed th e steady state solutions of MHD with time translationally invariant 
stochastic force terms. As shown in Sec. VIII, information about cascade directions also can be deduced from this 
framework. To our knowledge, this is the first time the RG has been used to predict cascade directions in turbulence. 
In that section, we examined directions for energy cascade in both hydrodynamics and MHD and helicity cascade in 
MHD. In the stability region of the trivial fixed point, in all cases and in any d the one loop RG analysis showed 
no cascade at asymptotic scales. On the other hand, in the stability region of the nontrivial fixed point, the RG 
analysis predicts for the energy in all cases and for all d, direct or normal cascades and for helicity in MHD also a 
direct cascade. For energy cascades in d = 2, 3 for MHD and d — 3 for hydrodynamics, the results of the RG analysis 
found here are consistent with general expectations. However, for energy cascade in d = 2 hydrodynamics and helicity 
cascade in d = 2, 3 MHD, some general arguments (Kraichnan 1967; Frisch et al. 1975) and simulations (Lilly 1969; 
Frisch and Sulem 1984; Sommeria 1986; Pouquet at al. 1976; Pouquet and Patterson 1978; Pouquet 1978) suggest an 
inverse cascade. There are no rigorous proofs in any of these cases. Nevertheless, the cases where the RG analysis 
differs from conventional expectations deserves further discussion. Below we identify some possible sources for the 
discrepancy and suggest ways they might be resolved. However first, let us consider the random force terms. The 
class of random force functions we consider are reasonably general in terms of covering a broad range of scales, and 
we do not believe the special form of the force is the main cause for this discrepancy. These have been taken to be 
power laws in wavenumber, characterized by a free exponent y that can either be positive, zero or negative. This is 
enough to cover most cases of spatially correlated noise and white uncorrelated noise (if y = 0). 

Returning to the discrepancies, there are several possible reasons that come to mind for these disagreements between 
our RG analysis versus the various simulations and general arguments. We will first list them and and then address 
each in turn. The main possibilities are: (a) intrinsic limitations of perturbative RG applied to MHD; (b) limitations 
of the one-loop RG analysis which may be corrected by higher loop calculations; (c) the RG analysis is correct and 
the inverse cascades suggested by simulations are a metastable effect in the unstable region of the nontrivial fixed 
point. 

Concerning point (a), we remark that it has long been known and criticized that since the formal (bare and 
unrenormalized) perturbation expansion parameter Ao = 1 in the Navier-Stokes equation is equal to unity, the 
(primitive) perturbation series based upon it are of dubious utility (see, e.g. discussion in McComb 1995, and the 
references therein). Analogous primitive perturbation series constructed for MHD would not escape this criticism 
either. In response, we point out that in the loop expansion carried out here, we have defined and expanded our 
perturbation series in powers of the amplitude of the noise fluctuations (which could be large or small). This is 
because the noise amplitude serves to count the number of loops in the series (represented by a graphical expansion) . 
An expansion in Ao is not the same as an expansion in loops. The parameter Ao is only a strict bookkeeping device, 
or "flag" , allowing us to keep track of the nonlinear adjective mode-mixing terms. Note moreover, that Ao is a 
dimensionful parameter, and it makes no physical sense to expand in terms of it. Meaningful perturbation series must 
ultimately refer to dimensionless expansion parameters or couplings. The actual perturbative expansion parameter 
for MHD is not Ao, but rather, g±, which as we have calculated, is less than unity at both the fixed points. Thus, it 
is even possible to have reasonably large fluctuations and maintain a small expansion parameter, (see Eq. (0)). This 
coupling g± is dimensionless, and provides a bona-fide expansion parameter. So, we do not believe there to be any 
inherent internal inconsistency in making use of a perturbative loop calculation of hydrodynamics or MHD, provided 
we correctly recognize that the true (and small) expansion parameter is not Ao but rather Eq. ([i"2"|). 

Turning to point (b), this has to do with possible limitations of the one- loop analysis itself. The number of loops 
considered is intimately related to the magnitude of the stochastic fluctuations. If this number is small, then one 
reasonably expects a low order loop expansion should be able to capture the salient physics in the problem. However, 
if the fluctuations are large, one may have to push the calculation to higher loops. Once again however, the noise 
amplitude is not a dimensionless parameter, and what counts loops, properly speaking, is the dimensionless ratio 
written in Eq. (^). Note that large effective viscosity and/or magnetic resistivity can counteract the presence of 
large stochastic fluctuations. The number of terms or graphs in the loop expansion grows geometrically with the 
number of loops, so ample physical justification would be needed before taking the calculational leap. A tempting 
and quite interesting, though exceedingly nontrivial, proposal would be to try and go beyond the loop expansion 
by using nonperturbative methods to make predictions of cascade directions, following Wilson's ideas for the exact 
rcnormalization group (Wilson and Kogut 1974). In any case, neither of these options have yet been pursued. 

Finally, coming to the point raised in (c), it may very well be the case that the RG analysis is correct in all its 
physical predictions despite the one-loop order, and there is no actual discrepancy with the numerical results. The 
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key to the resolution is in recognizing that the numerical simulations are always limited by the finite duration times of 
the computer runs and the finite size of the lattice, whereas the RG analysis here refers always to the asymptotically 
largest spatial and temporal scales, which are evidently not accessible to present day computer simulations. So, we 
are comparing phenomena (cascade directions) associated to very different spatial and temporal regimes. In other 
words, the inverse cascade might not be an asymptotic phenomena, but is linked instead to intermediate scale physics, 
and this is evidently what the simulations are seeing. One possible explanation is that simulations are seeing a low 
frequency mode that temporarily shifts energy to larger length scales but eventually shifts the energy to smaller length 
scales. The RG analysis in this paper would ignore such modes. Let us elaborate further on this suggestion. Both 
inverse and direct cascades are possible, but different cascade directions would come to predominate at different time 
scales. The complete dynamical evolution may occur in two phases: an earlier inverse cascade phase followed by a later 
direct cascade phase. If true, this two-phase behavior would be an example of a "cross-over" phenomenon, familiar 
from critical dynamics. Of course, we do not see any evidence for cross-over from our one-loop RG calculation, and 
only a higher-order (or preferably, non-perturbative) calculation could begin to shed some light on this possibility (see, 
point (b) above). To get a further hold on the problem, one could either (1) run the simulations for ever increasing 
times and larger lattices, or (2) study the RG away from the vicinity of the fixed points. In the latter case, we are 
"backing away" from the leading asymptotic limit and the terms in the equations of motion denoted as "less relevant" 
in Eqs. ( p7|) and (|29| ) would now have to be included in the RG analysis. 
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APPENDIX A: DYNAMIC GENERATING FUNCTIONAL FOR STOCHASTIC MHD 

The purpose of this Appendix is to set up and apply a powerful and extremely useful dynamic functional integral 
technique that leads to an efficient, systematic and organized perturbative loop-expansion for the solutions of the 
fully nonlinear magnetohydrodynamic equations Eqs. (^J) and (|To|), when the latter are subject to random forces and 
fluctuations. The ability to cast stochastic differential equations in terms of functional integrals opens the door to 
the calculation of many physically important quantities including the so-called effective action and effective potential 
(Hochberg et al. 1999b). These quantities can be used to investigate the onset of pattern formation and pattern 
selection and can be used to identify symmetry breaking states of the system. In this paper, we will use functional 
methods to deduce the Feynman rules needed to construct a compact graphical representation of the perturbative 
solutions of the coupled nonlinear MHD equations. These rules permit one to identify a distinguished class of diagrams, 
the so-called one-particle irreducible (1PI) graphs. The renormalization of these diagrams to a given number of loops 
will be used to derive the differential renormalization group equations Eq. ( |io| ) , that govern the scale dependence of 
the parameters appearing in the MHD equations. Additional details regarding the application of functional and path 
integral methods to stochastic differential equations can be found in (Martin et al. 1973; DcDominicis and Peliti 1978; 
Rivers 1987; Zinn-Justin 1996; Hochberg et al. 1999b) and references therein. In the following, vector quantities will 
generally be denoted in boldface (or sometimes with an arrow) and matrix operators indicated by means of the dyad 
(double-arrow) notation. 

Let QsoimPsoin denote a general solution of the coupled stochastic MHD equations Eqs. @ and ([To|). We pass 
immediately to a functional or path-integral representation for an arbitrary function F of this solution-pair by means 
of the trivial identity 

/ [VP][DQ] F(Q, P) d[P - P soln ] S[Q - Qsoin], (Al) 

where the path integrals are taken over all field configurations. The product of delta-functionals selects out the 
particular field configurations P = P so / n and Q = Q. so ;„ corresponding to the solutions of the MHD equations. They 
constrain the functional integrals which otherwise freely sum over the totality of all Q and P field configurations. 
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The solutions of the MHD equations depend of course on the random noise sources, so that P so in = P so i n [fjp : rfg], 
and Qsoin = QsolnlVPi Vq]- The solutions themselves can be regarded as functionals of the noise source. Ultimately, 
we are only interested in the solutions of MHD averaged over the random forces and random currents. This statistical 
average is denoted by angular brackets and is given by 

(F(Q,P)) = f [Vff\F(Q 8oln ,P 8oln )V[ff\, (A2) 

where the noise distribution functional V is assumed to be given (i.e. it is a Gaussian functional for Gaussian noise) 
and is normalized to unity 

J [DrJ P }[Drf Q }V{rjp,ri Q } = J [Vfj\P[fj\ = 1. (A3) 

We abbreviate the 2<i-component noise vector field by ff = (ffp,f}Q), where ffp and f]Q each have (i-components in 
a (i-dimensional space. The goal and objective of the following sequence of manipulations consists in choosing a 
particular form for F and explicitly performing all the functional integrals encountered along the way. The resulting 
functional then is used to generate all the statistically averaged products of solutions of the nonlinear MHD equations 
to any desired order in the expansion parameter Ao- The Feynman rules (see Fig. 1) can be read off directly from 
this final functional. 



We need to carry out a change of variables in Eq. (Al). To see how this works, consider first a simple example. 



Let f(x) = (fi(x,y), f2(x,y)) be a two-component vector- valued function of a vector variable x = (x,y) and suppose 
this has a unique root at xq, f(xo) = 0. Now make use of the functional generalization of the identity 

S(f(x)) = 6(h(x, y ))6(f 2 (x,y)) = /^~5°L (A4) 

det IlM ' y} 



d(x,y) 



where 

dfi,2(x,y) 
d(x,y) 

and 



dfi(x,y) 8f 1 (x,y) 

dx dy 

9f2(x,y) df 2 (x,y) 

dx dy 



(A5) 



S(x - x ) = S(x - x )S(y - y ), (A6) 



to express the product of delta functionals in Eq. (Al) as 



6[P - Psoin] S[Q - Qsoin] = J5[d t Pj + A P,„^ (Q t P n ) - 7+V 2 P, - 7 -V 2 Q, - (r) P ) 0l 

2/n w2 , 



x S[d t Qj + XoP^AiPtQn) - y+V A Qj - 7-V 2 P, - ( VQ ) 3 }. (A7) 
Here, the Jacobian functional determinant J is the determinant of the following 2 x 2-block array: 



(d t - 7+V 2 )<S jfe + X P jk di(Qi) -i-V 2 6 3k + X P jn d k (P n ) 
- 7 _V 2 ^ fc + X P jn d k {Q n ) (d t ~ 7+V 2 )<5 jfe + A P, fe 9z(P/; 



5(2-a?)5(t-t'). (A8) 



When the identity Eq. (A7) is substituted into Eqs. (Al) and (A2), the stochastic MHD equations Eqs. ( fl9| ) and 



( pp| ) appear explicitly under the functional integrals as constraints, and we see that the stochastic dynamics is built 
in from the outset. The next step is to replace each of delta- functionals appearing on the left-hand side of Eq. ( A7) 
with its (functional) Fourier integral representation. That is, for any (vector-valued) space-time field <&(x,t), 

$[*] = J pa] exp (i J d d xdt a- *), (A9) 

where a is the Fourier conjugate field, and has the same number of components as Also, in the above, the dot 
product indicates a sum over the repeated indices, so that <? ■ <& = crj&j, for j — 1, 2, • • • , d. Replacing each delta- 
functional by its corresponding Fourier integral introduces two conjugate field variables into the stochastic average Eq. 
(A2). Since there are two Elsasser fields, namely P and Q, each requiring a separate conjugate field, we distinguish 
each individual conjugate field by a superscript: a = (a p , a®). 
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At this stage what results is a multiple functional integral representation for the stochastic average of an arbitrary 
function F, which depends on the exact solutions of the original pair of MHD equations. That is 

(F(Q,P)) = J [VfJ\ [VP] [VQ] [Va p ] [Va Q ] F(Q, P)J V[rf\ 

x exp (i J d d xdt afldtPj + \ P jn di(QiPn) - 7+V 2 P, - 7-V 2 Q, - (r/ P ),]) 

x exp (i J d d xdt af[d t Qj + X P jn di(PiQ n ) - 7 +V 2 Q, - - (tjq)^). (A10) 

Following the standard practice in both quantum and statistical field theory, now choose a convenient form for F. 
In particular, let J p , , 3" P , 3 <jQ denote a set of four arbitrary vector source functions, one source for each of the 
Elsasser and conjugate fields, respectively. Then let 



F = exp (P • J p + Q • J Q 

and define the generating functional Z by 

Z[3 P ,3 Q ,r P ,r Q ] = (exp(P-J p - 



a J 



Q • J Q 



J" 



r Q )). 



(All) 



(A12) 



We have suppressed the omnipresent integrations over space and time to avoid clutter (they can be restored at any 
time). The source functions allow one to factor the non- linear (in this case, cubic) interaction terms out from under 
the functional integral, where they now appear as functional derivative operators acting on a simpler and exactly 
calculable generating functional. Using the functional generalization of the identity 



du g(u) exp(wJ) = ff(Jj) / du exp(uJ), 



(A13) 



it then follows from Eqs. (A1C), (All) and (A12), that 



Z[J p ,J Q ,J ffP ,J CTQ ] =cxp (i J d d 



, , & x, n r 8 8 

Xdt Xo - „ P PjnOl 



SJf 



x exp (i J d d xdt\ j^ T P jn di 



6J?SJP- 
S 6 ■ 
WW 



(A14) 



where Zq is the functional Z evaluated at zero coupling: Ao = 0. The functional derivatives have the property that, 
e.g., 



8Jk(x,t) 
6Ji{x',t>) 



S kl S(x-x')S(t~t'). 



(A15) 



Before continuing, special mention should be made regarding the treatment of the Jacobian functional determinant 
J that appears in Z. The determinant can be exponentiated by means of the identity detM = exp(trlnM), and 
it can be shown that the argument of the exponential is proportional to the formally divergent object 5 d (0). This 
happens to be an ultraviolet or short-distance divergence, and when renormalized by dimensional continuation, yields 
zero, S d (0) = (Zinn- Justin 1996). Thus, in this regularization scheme, the determinant renormalizes to yield a finite 
and field independent constant (unity), and so contributes nothing to the normalized correlation functions calculated 
from Z. We henceforth drop these determinants for the remainder of this discussion (by setting J — > 1). If one 
wants to investigate the short-distance or ultraviolet limit of turbulent MHD however, then it is preferable to keep 
these Jacobian determinants i ntact . For example, this is important when calculating the effective action and effective 
potential associated with Eq. ( A14 ). The retention of these determinants can be done with ghost fields, as described 
at length in Hochberg et al. (1999b), where a number of details pertaining to the calculation of Jacobian determinants 
can be found. To our knowledge, the calculation of J for MHD is an open problem. 

Up to this stage, we have assumed a completely general noise distribution functional V . For the remainder of 
this section, and throughout this paper, only Gaussian noise is considered. For Gaussian noise, Zq can be computed 
exactly and in closed form. In fact, we consider translationally invariant Gaussian noise, that is 



V[ff[ = A/" exp ( - i I d d xdt I d d x" dt 1 ? ?i (x, t) TT 1 (x -x',t- t') r,j (x' , t')) , 



(A16) 
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where M is a normalization constant and r=T is the symmetric covariance matrix, which will be specified below. By 
substituting this form of V into Zq, and completing the square in the conjugate field variable, which is a 2<i-component 
vector field a — (a p , o®), exact integration can be done immediately over the noise field^. 
Following this procedure, completing the square yields 

1 w-l 1 «-» «->-l <-> 1 <-> 

- -ff- r -ff- ia-ff = --[ff+i T -o] T - T -[ff+iT-a] — -a- T ■<?, (A17) 
where T denotes transpose. Next, the functional integration over the noise leads to 

Z [J,J CT ]= y [£><?] [£>#]cxp(- -ff. ? -ff + iff- G ■* + *■ J + a- J ff ) , (A18) 



where we have used the fact that the noise distribution function V is normalized to unity Eq. (A3). We have 
streamlined the notation even further by defining two-component fields for the Elsasser fields, their sources and their 
conjugate fields respectively as 4? = (P, Q), J = (J p , J^), and J <T = (J CT , J CT ). The (free) inverse response function 
in the space and time domain is defined to be 



In Eq. (A18) we have a quadratic form in the conjugate field, so we can integrate over this field exactly by 



completing the square as 



1 <-» 

- -ff- T -ff + iff- G -* + <r-J c 



-\{S-iY 1 -[Gq 1 ■&-iJ a ]) T - r ■(?-!? 1 -[Go 1 ■* -*J CT ]) 

- i(Go * * - U ff ) T - F 1 -(So 1 * - iJ CT ). (A20) 



Integrating over the conjugate fields yields (up to an overall and unimportant multiplicative constant) 
Z [J,J CT ]oc Ax>*]exp(- -(Go* -*-?:J ff ) T - r 1 -(Gq 1 •*- i J CT ) + *- j). 



(A21) 



Finally, by means of the identity 



(g 1 * iJ CT ) T - r 1 -(G 1 * - «J") - 

(* - i Go -J CT ) T • (Go' • ST* • Gq" 1 ) • (* - * Go -J"), (A22) 

we are led to make the constant shift change-of-variablc <fr' = <I> i Go •J' 7 , integrate over the remaining field and 
note that = [2?<I?]. Thus, up to an overall irrelevant constant prefactor, the free (Ao = 0) generating functional 

is given by 

Z [J P , J g , 3°", r Q ] = cxp Go -j) exp (ij • [Go ■ ? • Go] ■ j). (A23) 
This can be expressed in Fourier space as 

Z [J P ,JQ , J- P ,J^]=exp(z J ^^J-(fc,c).5 (-k,-cj)-3(-k,-c 

/ 1 f d^k dui 1 » -» 1 > — ' ' \ 

x exp y (2 ^ d+1 J(fc, w). Go (fc, w)- T (-fc, -u)> Go (-fc, -w) • J(-fc, -w)J . (A24) 



2 Gaussian functional integration is treated in any number of excellent texts devoted to (quantum) field theory. We recommend 
the excellent coverage given in (Ramond 1981; Rivers 1987; Zinn- Justin 1996). 
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From the free generating functional, it is straightforward to read off the bare response and bare correlation functions 
(or by taking the appropriate functional deri vatives of Zq with respect to the sources). In fact, the whole point in 



defining Z the way we did back in Eq. (A12) was so we could compute the statistical average of products of fields 
and conjugate fields. As we see, functionally differentiating Z with respect to any combination of sources brings 
down the corresponding factors of fields under the functional integral. Then, setting all the sources to zero after 
differentiating yields the correlation and response functions associated with this stochastic field theory. Let us now 
specify the explicit form for the noise or covariance matrix as 

[? (fc, L0)]mn = Do(k) P m „(fc), (A25) 

which is the product of a 2 x 2-block noise matrix Do(k) times the transverse projection operator. The presence 
of this projector enforces the solcnoidal character of the noise, and appears as well in the original noise correlation 
functions above Eqs. (^) and (CIC). The noise matrix can be used to model an arbitrary spatially correlated noise. 



To handle both spatial and temporal correlations in the noise, one would write Dc\( k, lu). In this paper, we consider 



spatially correlated noise as specified by 2 x 2-block D${k) given below in Eq. ( A2S ). The e ffects of random helicity 
(both kinetic and magnetic) can be studied by an appropriate and simple extension of Eq. (A25); this is introduced 
and analyzed in the main text. 

Thus, the bare response function in Fourier space is computed from 

(2 7r ) 2 ( d + 1 ) S 2 
(cr m (pi,wi)$„(p 2 ,W2))o = — ^r™ — xt<t(~ u 7 (-* dj-=J=o^o[J CT , J] 

= i(2n) d+1 5 d {p 1 +p 2 )S(uJi +w 2 )[Go (pa,wa)W (A26) 
The bare response function obtained by calculating the inverse ofEq. ( |A19|) is 

[Go {k,uj)] mn = Go(k,u))S mn , 

r it \ — ^ ( —iuj + -f+k 2 — 7_fc 2 A , 

y ° { H^ + fc 2 (7 + +7-)][-^ + fc 2 (7+-7-)] I -7-fc 2 -^ + 7+ fc 2 J' (A27 > 

where the convention is that momentum and frequency are taken to flow from the conjugate field to the physical field 
(otherwise, one must flip the sign of &). Next, the bare correlation function (in Fourier space) is computed from 

(27r) 2(d+1 ) 6 2 
= {2n) d+1 8 d (p 1 +p a )8(u 1 +u 2 )[Co (fi,«a)U, (A28) 

where, 

[Co (k,uj)] mn = C (k,uj)P mn (k), 

C (k, uj) = Goik, w)D (k)Go(-k, -to), 

D »^=\2k-VB 2k-vA)- (A29) 
To complete the specification of the full generating functional Z, we must also express the non-li near i nteraction 



verte x in Fourier space. This is obtained by Fourier-transforming the cubic terms appearing in Eq. (A10) or in Eq. 
( |A14| ) . Note there are two kinds of cubic interaction, one corresponding to the er p — P Q vertex and the other 
corresponding to the er Q — P — Q vertex (the distinguishing feature is the conjugate field). Of course, there is only 
one (bare) cubic coupling parameter Ao (which we set to unity at the end of the calculation) . This means corrections 
only to one of the two vertices needs to be considered, even though both bare vertices are needed in building up 
loop-expansions. Taking the convention that all momentum and frequency flow inward to the vertex, the bare vertex 
corresponding to the cubic interaction o-?(kx,u)\) — P n {k2,^2) — Qi(k3,uJs) as well as the other cubic interaction 
crf(ki,u-i) - Q n (k 2 ,u} 2 ) - Pi(k 3l uj 3 ) is 

i\ Pl, jn (ki) = t\ (fci),P jn (fci) - i A (fei)A-„ - fa)ifa) n \ (A3Q) 
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The vertex is independent of the external frequencies. Also, here it has been written so that it depends only on the 
momentum carried by the conjugate field(s), since momentum conservation implies k\ + fca + &3 — Cl- 
in summary, in this Appendix the stochastic PDEs of MHD Eqs. ( |l9| ) and (20) have been e xpres sed in terms 
of a functional integ ral representation. From this representation, the bare response function Eq. (A27), correlation 
functions Eq. (A2£), and the bare non-linear interaction vertex Eq. ( A30| ) have been identified and calculated. As 
explained earlier, this set of functions constitute the basic elements (see Fig. 1) upon which a diagrammatic expansion 
(see Fig. 2) of the generating functional can be systematically developed. 



APPENDIX B: RESPONSE FUNCTION RENORMALIZATION 



This Appendix carries out the one-loop renormalization of the response function. Transcribing the diagrammatic 
equation in Fig. (2a) for the one-loop corrected response function into the corresponding analytic expression by means 
of the Feynman rules in Fig. 1 yields 



[G {k,Lo)] mn = [Go [k,u))} mn -5Aq[Go (k, bj)] mm [P;, fis (k) + P s ,ni{k)]Irni s (k, to)[Go (k,uj)]r 
where the symmetry factor for the one loop diagram 5 = 4 and the loop-integral is given by 
*> d d q f°° dfl 



<(k,u) 



(2ir) d 7-oc 2vr 



\Pm,ir{k - q) + P r ,irn(k - q)]Qo(k -q,LJ- Q)C (q,fl)P rs (q). 



(Bl) 



(B2) 



The precise meaning of the momentum-shell integration will be given below. Using th e fact orization of the response 
function into the product of a block matrix times a Kronecker delta as given in Eqs. ( A27 ) and (Bl), the following 
equation is obtained for the one-loop correction to the inverse response function, 



Q- l {k,ui)8 mn = G^\k,Lu)S mn + SXl[P ltns (k) + P s ^ nl (k)]I mls (k,u;) + 0{4), 
where the renormalized 2x2 inverse block is parametrized as 



g-^k^) = 



7 5fc 2 



l<k 2 
-ilu + 7+ k 2 



(B3) 



(B4) 



This structure follows from Fourier transforming Eq. ( A19). The infrared (i.e., large-distance and long-time) renormal- 
izability of the response function implies that the renormalized or corrected response function Eq. ( B4 ) must have the 
identical mathematical structure as its bare counterpart Qq ■ It therefore has the same frequency and momentum de- 
pendence as the bare response function and contains the same number of parameters. Here 7 ± = 7^ (A, B, 7+ , 7_ , Ao) 
are computable functions of the indicated bare parameters which will be determined below. 

Thus the renormalization of the response function yields the individual renormalization of the two independent 
viscosities 7+ , 7_ . The renormalization of the two viscosities requires some amount of careful calculation. We will go 
through the major points step by step. To proceed we need to explicitly carry out the internal frequency integration 
that appears in the loop integral Eq. (p32]). From Eq. (B4), we see that the viscosities 7± are the coefficients of the 
k 2 terms, so we can set the external frequency uj to zero from the outset. This frequency integral is straightforwardly 
calculated by means of the calculus of residues. We close the contour in the lower half plane, where there are two 
simple poles, and obtain 



dVL 
2^ 



g ( k -q,-n)c (q,n) 



Do(q) 



q ■ k 

-I- 



{} + q -f + o { k 2 

q ■ k 





1 -l" 


1 


'l 1 


){(7+~7-) 2 


-1 1 


+ (7+ + 7-) 2 


1 1 



A — B 
(7+ -7-) 2 



-1 1 



A + B 

(7++7-) 2 



1 1 
1 1 



(B5) 



The result displayed in the first line holds for all spatially correlated Gaussian noise D (q). We should point out this 
was evaluated assuming a positive viscosity, v > 0, and a positive resistivity, vb > 0. A change in sign will change the 
location of the poles and correspondingly modify the frequency integ ratio n. The second line results from substituting 
the particular form for the noise spectrum used in this paper, Eq. (|A29|). If one wants to treat noise with temporal 
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correlations, i.e., noise spectra of the form Do(q, SI), then the frequency integral must be evaluated in a case-by-case 
fashion, and will generally involve both p oles and branch cuts. 

Substituting this result bac k in to Eq. (JB2|), we turn next to the remaining and final integration over the internal 
loop momentum. From Eq. (B4) the renormalization of the parameters 7± requires that we expand the one-loop 
contribution up to second order in the external momentum. The projection operator in Eq. ( p33| ) is linear in fc, 
so the loop integral only need be expanded to first order in external momentum. It is important to note that this 
integral depends on the external momentum not only via the integrand but also through its domain of integration. 
This is because, when integrating over the momentum shell, we must ensure that all combinations of internal and 
external momenta appearing in Eq. (B2) remain within the band or shell of momenta to be integrated over. In the 
case of the present integral, this means we must integrate over the intersection of the domains A/s < \q\ < A and 
A/s < \k — q\ < A, where s > 1 is a measure of the fraction of modes to be eliminated in the coarse-graining. To first 
order in k, the second inequality can be written as A/s + k cos 9 < q < A+ k cos 9, where 9 is the angle between k and q. 
There are two cases to consider: when cos 6* > 0, the intersection of the two intervals can be expressed as the difference 
[A/s, A] — [A/s, A/s + k cos 9), and when cos 9 < 0, the intersection can be written as [A/s, A] — [A + k cos 9, A]. This 
means the complete momentum- shell integration valid up to 0{k 2 ) can be written as 

idy r . r A r A/s+kcos9 r A \daa d - 1 

' - f dn d - - m_ + 0(fc2) . ( B6) 



(2ir) d J Wa/s J a/s JA+kcose' ( 27r ) 

Furthermore, note it is convenient to carry out the integrations over the momentum modulus before working out the 
angular integrations over the element of solid angle dVL d . 

From Eqs. ( jB5| ) and ( ]B2] ), the remaining integration to be carried out after the frequency integration is performed 

is 

1= f ^[V^k-^+VrM^-qWrsi^q-V-^l + ^ + Oik 2 )), (B7) 
= h + I 2 + h, (B8) 
where we break the full integral / down into somewhat smaller and more manageable pieces as 
"> d d q 

W) 



h=j j^- d q-y- A (-q m V ls (q))+0(k*), (B9) 

^=/ > (0^ 4 (-^®)^+°( fc2 )> ( bi °) 

^q~V-^k m P^-(k s -k-qp^ + k r P lr ,MPr S (q))+0(k 2 ). (Bll) 



h = 



The following projection operator product expansions proved to be useful in arriving at these expressions: 
Pa.Uk q)P dc (q) = (k a - q a )[6 bd - - ^ [k d - ^-^) + 0(k 2 ), 



q-k_qd s 

aM {k - q)P ad {q) = k a P bc (q)Pad(q) + 0(k 2 ). (B12) 



The integrand of I\ is independent of k and the contribution from the interval [A/s, A] vanishes upon angular averaging. 
Only the latter two contributions in Eq. QB6| ) will contribute to the angular averaging and we find that 

I, - [(^r- 4 - A— 1 {l^y d k m S ls - \ d[d+ % {2n)d ik m S ls + k l5ms + k M ) , (B13) 



where S d is defined in Eq. (E2). The integrand in I2 is already linear in fc, so it suffices to work out the contribution 
from the interval [A/s, A]. This yields 

h = (a 2 T\ h - ( B14 ) 
{d-y-4) 

Lastly, the integrand in ^3 is also linear in k. Integrating over [A/s, A] and averaging over the angles yields 
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h = 



[A d-y-4_ { A )d -y-4 ] {d _ 1)Sd 



(d-y-4) 



d{2ir) d 



-k m Si s + 



(d-3)S 
d(2ir) d 
2S d 
d(d + 2){2n) 



ks^lm 

-^[{k m Si s + kiS ms + k s 5 m i)j . 



(B15) 



Various identities needed for the angular integrations are collected i n Eq . {E2V T his complete s the c alcul ation of the 
one-loop integral. Combining all the results from Eqs. (B7), (p^), ( B10| ), ( P 1 1 ) , (B13) ,(B14), and (B15), we obtain 
the one-loop matrix equation 



Q- x {k,Q)8 mn = gQ\k,0)S mn + S 



4 d(d + 2){d-y-4){2n) dl V J 



f A — B 


1 -l" 


A + B 


1 1 


l(7+-7-) 2 


-1 1 


+ (7++7-) 2 


1 1 



x k P mn {k). 



(B16) 



We can project out individual equations for each of the two viscosities by multiplying Eq. (B16) through by P m i(k) 
and using the idempotency of this projection operator, P mn (k)P m i(k) = P n ;(fc). Thus we find that each viscosity 
renormalizes as 



7< = 7+ +s^-\A d - y - 4 - (-) d -y- 4 ]. 

I+ 7+ 4 [ [ s> 1 d{d + 2)(d - y - 4)(2ir) d I (7+ - 7-) 2 (7++7-) 2 



(d 2 - y - 4)S d 



(A-B) 



(A + B) 



A 



A. 



= 7 _ +l s^[A d -'«- 4 - (_) d -a-4 



(rf 2 -y-4)5 d 



-(A -5) (A + B) 



d(d + 2)(d-y-4)(27r) d l( 7+ -7-) 2 (7+ +7-) 



(B17) 
(B18) 



Individual equations for the sum (7+ + 7- ) and difference (7+ — 7 _ ) can be easily obtained by taking the sum and 
difference, respectively, of Eqs. (B17) and (B18). This yields the individual renormalizations of the fluid viscosity v 
and magnetic resistivity vb- 



APPENDIX C: NOISE SPECTRUM RENORMALIZATION 

From the diagrammatic representation of the one- loop correction to the noise spectral function in Fig. (2b), (the 
noise spectral function is obtained by amputating the external legs from the correlation function) which is a 2 x 2 block 
matrix, we obtain the following expression for the one-loop noise spectrum in terms of the bare response, correlation 
and vertex functions, 

D{k)P rs {k) = D (k)P rs (k) +SA 2 [P ; , r „(fc) +P n ,rl$)]Il m ,nj(k,u)[Pm,sj(k) +Pj,«n(% (CI) 

where the symmetry factor in this case is S = 2 and the associated one-loop integral is given by 

f y d d q f°° dfl <-> <-> 
Ii m ,nj(k, w) = J j^-^ J — [Co (k - q, u - n)]; m [Co (q, to)] n j 

> d d q f°° dtt 



f drq f°° dVt - - 
J (2ttP J ^ P irn(k-q)Pn 3 (q)C (k-q,uj-n)Co(q,tt), 



C (k — q,w- Q,)C (q, Q) = Go{k - q,u - fl)D (k - q)Qo(q- k,Cl-cj) 

xg Q (q,n)D (q)g (-q,-n). (C2) 

The infrared renormalizability of the noise spectral function implies that the r cnorm alized spectral function must have 



the identical mathematical structure as its bare counterpart Do(k) in Eq. ( A29). It must therefore have the same 



frequency and momentum dependence as the bare noise spectral function and contain the same number of parameters, 

n( ,._(2k-*A< 2k-vB<\ f „. 

u ^>~\2k-yB< 2k-yA<J- 

Here A K = A < (A, B, 7 + , 7_, Ao) and B K — B < (A : B, 7+, 7_, Ao) are calculable in terms of the bare parameters and 
will be determined below. 

We are interested in evaluating this loop integral in the hydrodynamic limit (i.e., for small k and small lo). However, 
there is already an overall quadratic factor of k 2 multiplying the loop integral arising from the product of individual 



3G 



projection operator factors appearing in Eq. (CI). As can be seen from Taylor expanding about k = 0, the loop 
integral itself can only contribute a constant plus positive integer powers of external momentum k in the hydrodynamic 
limit. Thus, the complete one-loop correction in the long-wavelength regime is a power series in k, and it starts off at 
quadratic order. For the type of noise of interest to us, we now show that this one loop correction is irrelevant in the 
hydrodynamic regime. In fact, by setting k and ui to zero inside the loop-integral Eq. ( |C2| ), it gives 

r> d d a r°° dn 

Wi(0,0)= / — / —P lm (-q)P nj (q)C (-q,-n)C (q,n) 
J (2ir) a J_ oc 2n 

° dn 

— C (-q,-Q)C (q,n) 



f 


q d ~ l dq 


/A/a 


(2n) d 


f 


q^dq 


JA/s 


(2n) d 



dn d P lm (-q)P nj (q) / Z—JLM(q). (C4) 



Here 

\A(r,\ — I 

M 2 (q) Mi{q) 

is the 2x2 block- matrix fu ncti on of the modulus q = \q\, which results from integrating the product of the bare 



Ml,)-(m iJM'l, (C5) 



correlation functions in Eq. (C2) over the internal loop frequency. This can be straightforwardly calculated by means 
of residues (there are four simple poles in both the upper and lower half-planes, respectively), but we need not do so 
for the purposes of the present discussion. Note that the integration over the sphere factors out from the modulus and 
frequency integrations. This is because the product of projection operators depends only on angles, while the product 
of the correlation block-matrices depends onl y on frequency and momentum modulus. Integrating the product of the 
projection operators that appear in Eq. ( p4| ) over the unit sphere using the identities in Eq. (|E^) yields the result 

J dQ d Pi m (-q)P nj (q) = ^ {{d 2 - 3)5i m 5 nj + [Si n 6 m j +%<W])- (C6) 

After using the properties of the projection operators to further simplify the final expression, it leads to 
(k)+P n ,rl(k)]Ii (0,0) [P (k)+P jtSm (k)} = 

The one-loop correction is of ord er Q (k 2 ) and is proportional to the same projector as the bare noise spectrum. 



Inserting the above line into Eq. (CI) gives 



D(k) = D (k) + S\% 2k' x ( ^ff 2)( £)d ) j As ^ dqM ^ + °( fc3 )- ( C8 ) 



Referring to the form of the noise spectrum matrix in Eq. ( A29 ), we see that the one-loop correction to the noise 
spectrum is irrelevant in the long- wavelength limit provided y > — 2 (The case y = — 2 is akin to Model A of Forster, 
Nelson and Stephen 1977). The Taylor expansion of h m ,nj{k, 0) about k — (taking into account contributions from 



the domain of integration as well, see discussion surrounding Eq. (B6) ) will just introduce positive powers of the 



external momentum. Thus we conclude that at one- loop as k — > and for y > — 2 

D{k) = D (k) 

=> A< = A, (C9) 
and B< = B. (CIO) 

Note that for noise spectra having y < — 2, the one-loop correction will induce relevant operators in the hydrody- 
namic limit. For example, if y — —3, the one-loop correction goes as 0(k 2 ) and in fact dominates over the the initial 
bare spectrum in the infrared. To be able to renormalize consistently in these other situations requires adding the 
appropriate dominant relevant operators to the noise spectrum. For the case y = —2, one needs only complete the 



calculation sketched in Eq. (|C4|) and insert this into Eq. (C8). In this case, one would find that both A and B receive 
nonzero corrections in the hydrodynamic limit. 
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APPENDIX D: VERTEX RENORMALIZATION 



The diagrammatic expansion for the one-loop correction to the vertex function, or coupling constant, is depicted in 
Fig. (2c). As can be seen from this figure, there are three distinct one-loop diagrams contributing to the correction. 
Momentum and frequency conservation implies that the vertex function can depend at most on two independent 
external momenta and two independent external frequencies. With the momentum and frequency flow assignment 
taken as shown in Fig. (2c), where momentum and frequency enter via the conjugate field and exit via the two 
physical fields, the equation for the one- loop vertex correction is 



A < P;,„ m (fc) = X Pl,nm(k) 



5Ao[Pj,nt(fc) + Pt,nj(k)]Iimjt{k,k 2 , k - k 2 ;uj,uj 2 , u - w 2 ), 



(Dl) 



where the diagrammatic symmetry factor is S — 4. Here A < = X < (A, B, 7+, 7_, Ao) is the renormalized coupling 
parameter, which in principle can depend on all the bare parameters appearing in the equations of motion. The 
loop-integral representing the sum of the three types of (amputated) triangle graphs can be written as 

Ilmjt ( k, k 2 , k - fc 2 ; W, U>2, tO-U>2) = 
r> d d q f°° dQ r 

try J2^y J ^ \ c o(<f, H)Pjk(q)\Pk,ri{k2 - q) + PLrk(k 2 - q)]Go(q- k 2 ,Q - Lu 2 )8 rs 

x [P s , P m(fc - q) + ~P m , P s(k - q)]Go(q ~ «, f2 - w)8 tp 
+ Qo{-q, -Q)Sik.\PLkr(q) + P r .ki(a)]Cc)(a — fc 2 , fl - LJt)P rs (a— fc 2 

+ g (-g, -n)s Jk [p l:kr (q) + p r ,«(§)]g (*2 

x C (k-q,u-tl)P tp (k-q)}, (D2) 

which follows from translating the diagrams into their associated mathematical expressions using the basic elements 
derived above for the bare response Eq. (A27), bare correlation Eq. ( A29Q and bare vertex functions Eq. ( A3C| ). 
The trace (tr) is taken over the product of the 2x2 block matrices. We are interested in evaluating this integral 
in the hydrodynamic limit (i.e., for small k and small u>). As there is already a factor linear in k multiplying the 
loop integral, and the bare vertex is itself linear in momentum, we can immediately take the limit k — fc 2 — * and 
uj = ll>2 — > inside the loop integral right from the outset and focus attention on li m jt(0, 0, 0; 0, 0, 0). This results in 
a tremendous simplification. The properties of the projection operators allow one to further reduce products such as 



-q, -tt)5j k [Pi,kr(q) + Pr,ki(o)] c o(q- fca,fJ- w 2 )P rs (g- 

s, pm {k -q) + P m:PS (k - q)]Ga{q- k,Q- u))8 tp 



q,u>2~ tt)S rs [P P:Sm (q - k 2 ) 



,(<?- k 2 )\ 



[Pi,jr(q) + Pr,ji(q)][Pp,rm(q) + Pm,rp(q)]Pt P (q) = qiq m Pjt(q) = q hih m Pj t (q), 

and likewise for the other products that one encounters in Eq. (|)^). 
The angular integrations immediately can be computed by using 



dVL d hih m P jt (q) = 



with the result that 

A < P( ! „„ l (fc) = X Pi, nm (k) + Ajj( 



d(d + 2) 



(d + l)5j t 8l m - [Si 3 S mt + SitSjmij ■ 



(D3) 



(D4) 



-25, 



tr 



d{d + 2) 

l+1 d q r°° dn 

2^ 



l.nrn 



(k) 



(2n) d 

g (-q, -n)g Q (-q, -o)c (-«r, -n)} . 



Pm,nl(^)] X 

{Co(q, n)G Q (q, n)G Q (q, SI) - Go(-q, -O)C (g, n)Go(q, O) 



(D5) 



The frequency integral can be performed by the method of residues. The integrand has both simple and double 
poles in the complex frequency plane; the integration contours can be closed in either the upper or lower half-plane. 
The first and third factors appearing under the frequency integral have simple poles in the upper and lower half-plane, 
respectively. These separate contributions can be neatly combined to yield 



dQ 
2^ 



(C (q, n)Go(q, n)Go(q, Q) + G (-q, -n)G (-q, -Q)C (-<f, 



D (q) 



(7+ - 7-) 3 



-1 1 



(7+ +7-) S 



(D6) 
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On the other hand, the middle factor has double poles in both the upper and lower half-planes and its contribution 
to the frequency integral turns out to yield 



f°° dtt 

/ 7T - n)Co(£ n)g (q, ft) = 



Do(q) 



1 



(7+ - 7-) 3 



1 -1 
-1 1 



1 



(7+ + 7-) 2 



1 1 
1 1 



(D7) 



Subtracting Eq. (D7) from Eq. ( |D6| ) as per Eq. (D5) yields a null result for the frequency integral that holds 
for all spatially correlated Gaussian noise Do(q). Thus the vertex coupling parameter does not renormalize in the 
long- wavelength limit and we have that 



A . 



(D8) 



This non-renormalization of the nonlinear coupling parameter holds also for the Navier-Stokes equation driven by 
arbitrary spatially correlated Gaussian noise and is ascribed to the Galilean invariance of the stochastic NS equation 
(Forster et al. 1977). The importance of Galilean invariance for the Burgers equation (and the KPZ) equation is 
further addressed in Medina et al. (1989). We suspect that the non-renormalization of the nonlinearity holds to all 
orders in the loop expansion for stochastic MHD. To establish this one could investigate the renormalization of the 
vertex coupling via suitable Ward identities, as was done for the case of the KPZ equation (Frey and Tauber 1994). 
We leave that problem for a separate discussion elsewhere. Temporally correlated noise breaks the Galilean invariance 
and can ( and usually does) lead to non-trivial renormalization of Aq. 



APPENDIX E: INTEGRATING OVER THE D-SPHERE 

The various angular integrations and angular averages that are needed in the renormalization group calculations 
are collected here. The area element of the unit sphere in d-dimensions is 

dfl d = dcf) sinflidfli sin 2 6 2 d6 2 ■ ■ • sin (d ~ 2) 6> rf 2 rf6» rf 2 , (El) 

where < <fi < 2ir and < 8 } < ir, for j = 1, 2, • • • (d — 2). Let Sd represent the area of the unit d-sphere, fij denote a 
unit vector in the j-th direction, and T denote the Gamma function. Then 

27T d / 2 



dfld TliTljTlfiTliyi — _|_ 2) (^^^ mn ~^ ^im^jn ~t~ ^in^jrn) • (E2) 

The angular average of the product of an odd number of unit vectors over the unit sphere vanishes identically. 
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FIG 1: Diagrammatic (Feynman-type) rules for the perturbation theory expansion of MHD. 
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FIG 2a: One-loop equation for the response function. 
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FIG 2b: One-loop equation for the noise spectrum. 





